Process and signature of initiation on a finite fault system: a spectral approach

Summary

We present a numerical approach, based on a spectral analysis, for the initiation of the unstable slip on a finite fault region. First we focus on one fault model. We study the relationship between the weakening parameter and the largest positive eigenvalue of the dynamic spectral problem. Since the numerical approach based on the integral equations proposed by Dascalu et al. (2000) is appropriate only for small eigenvalues we use a finite element method which permits accurate computations for large eigenvalues. We show the relation between fault length and the first eigenvalue that governs the duration of initiation duration. The value of the weakening rate can be evaluated from the strain field in the elastic medium over a domain of confidence. A specific pattern of deformation represents the signature of the initiation phase. The spectral analysis and the numerical methods used for the single fault model remain valid for more complex fault systems. The interaction between two faults is examined. Finally, we introduce the concept of spectral equivalence between a heterogeneous fault system and a homogeneous fault with renormalized friction law.

1 Introduction

Earthquake prediction is actually the Holy Grail of seismologists and the existence of precursory deformation is a big issue. Laboratory experiments on friction (Dieterich 1979; Ohnaka et al. 1987; Ohnaka & Kuwahara 1990) point out the existence of a phase of slow motion that precedes the propagation phase, leading to the rupture. They pointed out that this preseismic slip should be recognized as a manifestation of the initiation process preceding the dynamic rupture. Throughout this paper, we will refer to this phase as the rupture initiation or initiation. Ohnaka et al. (1987) showed through high resolution friction experiments, the relevance of the slip dependent friction law in the initiation process. They showed the existence of a characteristic length scale L_c, the critical slip. The shear stress degrades from the initial state τ_i to a residual dynamic level of friction τ_d with the ongoing slip. They also demonstrated the existence of a nucleation zone over which the nucleation process occurs. The length of this zone can be related to the constitutive parameters of the friction law.

On the basis of these experiments, Campillo & Ionescu (1997) studied the initiation phase (Ida 1972) of an unstable antiplane shear process on an infinite fault under linear slip-weakening friction. They gave an analytical expression of the slip, derived from an eigenvalue analysis. Considering only the part of the displacement associated with positive eigenvalues, they defined the dominant part w_d(t, x, y) as follows:

where α is the weakening rate on the fault (y = 0) divided by the rigidity, w₀(x, y) and w₁(x, y) the initial perturbation in displacement and velocity respectively. α is given by the following relation:

where G, L_c, μ_s, μ_d and S are respectively the rigidity modulus, the critical slip, the static friction coefficient, the dynamic friction coefficient and the normal stress. Eq. (1) shows that w_d(t, x, y) is characterized by an exponential growth with time and a simple exponential dependence of the slip distribution with respect to the coordinate perpendicular to the fault. In response to a small perturbation (w₀, w₁) the instability will develop in a limited spectral domain. The limiting wavenumber is a function of the friction law and the elastic properties. The dominant part gives an accurate description of initiation. This is illustrated by the numerical comparison performed by Campillo & Ionescu (1997). They demonstrated that the dominant part is indistinguishable from the complete solution, except in the very beginning of initiation, when the slip is of the order of the perturbation and therefore negligible with respect to the final amplitude of the solution. The expression of the dominant part associated with positive real eigenvalues is obviously noncausal, as remarked by Knopoff et al. (2000). Indeed, since all the propagative terms are omitted, the dominant part is not expected to be causal. During the initiation the dominant part grows exponentially with time, while the propagative terms scale with the perturbation. This means that the mismatch between the full solution and the dominant part at the causality limit also scales with the initial perturbation. In the context of earthquake study, we suspect the early triggering event to be small and in practice undetectable. For this reason violation of causality in eq. (1), that is not a mathematical problem, is not important from a physical point of view either, since the measurable slip evolution is almost perfectly described by the dominant part.

The slip during the initiation phase starts very slowly and gradually increases with time. The transition from the initiation phase to the propagation phase occurs at the point on the fault which has slipped by L_c, the critical slip. At this time the crack front begins to propagate along on the fault. The duration of the initiation phase is approximated by:

where l is the half width of the perturbation, b its distance to the fault, c the S wave velocity, W₀ and W₁ the weighted average of the initial perturbation (Campillo & Ionescu 1997). During initiation the slip grows on a zone whose characteristic length is given by:

Note that l_c is the only characteristic length scale associated with the infinite fault model. Dascalu et al. (2000) studied the initiation process on a 2-D finite fault of length 2a. They introduced a non-dimensional parameter β = a.α, and performed a stability analysis to determine a universal constant β₀ that depends only on the geometry of the antiplane problem. This constant was found as:

and defines quantitatively the limit between stable (β<β₀) and unstable (β>β₀) behaviour of the fault. For a given weakening rate on the fault, the minimum length for an unstable behaviour is:

The critical fault length 2a_c for unstable behaviour is a completely different concept from the critical crack length l_c, defined in the infinite fault context. As a consequence, when the fault length 2a is much larger than 2a_c, the pertinent characteristic length scale is l_c. The ‘free’ initiation process, as it would develop on an infinite fault, provides a good description of what really happens on the finite fault. The finite fault effect arises when the fault length 2a→2a_c. One can feel that this time, the ‘free’ initiation process (infinite problem) is unlikely to describe the real initiation process. The fault finiteness will have an effect on the initiation process. This effect was shown by Ionescu & Campillo (1999). The main effect is that the duration of initiation is greatly increased as the fault length 2a tends to 2a_c. Indeed, as it follows from Dascalu et al. (2000), the duration of initiation is T_c≈λ₀⁻¹, where λ₀² is the largest eigenvalue, and long initiation duration is expected when β→β₀. The goal is to determine how strongly the initiation process is affected by the fault finiteness.

Since the numerical approach based on the integral equations proposed by Dascalu et al. (2000) works only for small eigenvalues (i.e. β close to β₀) we present here a numerical approach to the initiation of unstable slip on a 2D finite fault, based on a finite element method (FEM). We compute λ₀, the largest positive eigenvalue of the dynamic spectral problem, as a function of β for the finite fault model. We study the asymptotic behaviour when β→∞ to deduce the infinite fault model. We compute the associated eigenfunction Φ₀(x, y) and the following ratio:

where ∂_y is the partial derivative with respect to the perpendicular coordinate. This ratio represents the information on the weakening parameter (Campillo et al. 2001). In the infinite case, this ratio exactly equals α everywhere in the elastic space. In the case of the finite fault zone, we point out the existence of an initiation pattern that develops in the elastic body. This pattern is associated with the unstable behaviour of the fault zone. Moreover, we show the existence of a domain over which γ is nearly constant. More precisely, in case of a single fault segment, this ratio equals the weakening rate α over a so called domain of confidence, whose size is an increasing function of α. Hence, the weakening parameter that is defined only on the fault can be estimated during the initiation phase outside of the fault, inside the bulk. This analysis remains valid for more complex fault systems. We first consider a system composed of two identical segments. We show that β₀ decreases as the distance between the two segments also decreases. This proves that fault segments may interact during the initiation process. Then, we consider a more complex fault system, composed of 11 identical fault segments, separated by rigid barriers. Each fault segment is characterized by the same weakening parameter β_local. We show the existence of the initiation pattern associated with each of the segments (locally, the fault system is reduced to only one fault segment). We also point out the existence of a global initiation pattern that is associated with the unstable behaviour of the whole fault system. We compute the dynamic spectral problem and derive the largest eigenvalue of the heterogeneous system. With the help of the relation between λ₀ and β computed in the one fault model, we get the corresponding equivalent weakening parameter β_equivalent. A simple homogeneous fault with a weakening rate β =β_equivalent will have the same largest eigenvalue λ₀, and therefore the same general average behaviour during the initiation.

The global initiation pattern associated with the unstable behaviour of the whole system is characterized by a wide domain over which the ratio −{[∂_yΦ₀(x, y)][Φ₀(x, y)]{ = β_equivalent/a. Locally, close to the fault, each individual fault segment has its own domain of confidence over which we can measure the local weakening rate. The interaction of all segments produces a global pattern of deformation at some distance from the fault system. This collective domain of confidence provides an estimation of β_equivalent. There is a complete spectral equivalence, in terms of the initiation process, between a complex fault system with β =β_local and a simple homogeneous single fault model with β = β_equivalent.

2 Spectral approach of rupture initiation on a finite fault

The method presented hereafter is basically the same for an infinite fault (Campillo & Ionescu 1997) or for a single finite fault (Dascalu et al. 2000). For convenience we give the principal steps of the spectral method in this section.

2.1 Problem statement

We consider the 2D antiplane shearing on a bounded fault region Γ_f included in the plane y = 0 in an homogeneous linear elastic space. The fault region can be composed of a set of simple faults on which the contact is described by a slip dependent friction law. We assume an antiplane configuration, that is the displacement field is 0 in directions Ox and Oy and that u_z does not depend on z. The displacement is therefore denoted simply by w(t, x, y). The elastic medium has the shear rigidity G, the density ρ and the shear velocity c = √G/ρ. The nonvanishing shear stress components are σ_zx = τ_x^∞ + G∂_xw(t, x, y) and σ_zy = τ_y^∞ + G∂_yw(t, x, y), and the normal stress on the fault plane is σ_yy = −S(S>0). The equation of motion is

for t>0 and (x, y) outside of the fault Γ_f. The boundary conditions on Γ_f are

if ∂_tδw(t, x)≠0 and

if ∂_tδw(t, x) = 0, where δw(t, x) = 1/2[w(t, x, 0+)−w(t, x, 0−)] is the half of the relative slip and μ(δw) is the coefficient of friction on the fault region. The initial conditions are denoted by w₀ and w₁, that is,

Since it is our intention to study the evolution of the elastic system near an unstable equilibrium position, we shall suppose that τ_y^∞ = Sμ_s, where μ_s = μ(x, 0) is the static value of the friction coefficient on the fault. We note that taking w as a constant satisfies eqs (8)–(11); hence w≡0 is a metastable equilibrium position, and w₀, w₁ may be considered as a small perturbation of the equilibrium. We shall assume that the friction law has the form of a piecewise linear function:

where μ_s and μ_d(μ_s>μ_d) are the static and dynamic friction coefficients, and L_c is the critical slip. Let us assume in the following that the slip δw and the slip rate ∂_tδw are non-negative. Bearing in mind that we deal with a fault plane and with the evolution of one initial pulse, we may (for symmetry reasons) put w(t, x, y)=−w(t, x, −y), hence we consider only one half-space y>0 in eqs (8) and (12). With these assumptions, eqs (9)–(11) become

where α is given by eq. (2). Since the initial perturbation (w₀, w₁) of the equilibrium state w≡0, is small, we have w(t, x, 0+)≤L_c for t∈[0, T_c] for all x, where T_c is a critical time for which the slip on the fault reaches the critical value L_c at least at one point, that is, sup_x∈Γfw(T_c, x, 0+)=L_c. Hence for a first period [0, T_c], called in what follows the initiation period, we deal with a linear initial and boundary value problem (8), (12) and (16).

2.2 The nondimensional problem and its spectral expansion

In order to obtain a non dimensional formulation let us introduce a the characteristic length i.e. we put

and we introduce the following nondimensional constant

Let us use the same notation Γ_f for the fault region in these new coordinates (i.e. x₁∈Γ_f if x∈Γ_f). From eqs (8), (12) and (16) we deduce

Let us consider the following eigenvalue problem connected to eqs (19)–(22): find Φ:R×R ₊ →R and λ² such that ∫^+∞_−∞ ∫^+∞₀ Φ²(x₁, x₂)dx₁dx₂ = 1 and

Since we deal with a symmetric operator we have real-valued eigenvalues λ², i.e. λ is real or purely imaginary. This symmetry property is specific to the slip dependent friction law used here.

The solution of eqs (19)–(22) can be generically written (in its spectral expansion) as:

where w^d is the ‘dominant part’ and w^w is the ‘wave part’, given by:

where (λ_n², Φ_n) are the associated eigenvalues and eigenfunctions of eqs (23)–(25), and

the projections of the initial data on the eigenfunctions. Let N be such that

The use of the expression of the dominant part w^d instead of w leads to a solution in which the perturbation is severely smoothed by the finite wavenumber integration. However, after a while, the propagative terms become negligible and the shape of the slip distribution is almost perfectly described by the dominant part.

2.3 Stability analysis; the static eigenproblem

One can easily remark that w≡0 is a stable position if λ₀² < 0 (i.e. N = 0). In this case the dominant part w^d vanishes and the system has a stable behaviour. Hence it is important to obtain a condition on β, that determines the positiveness of the eigenvalues λ². Since β is nondimensional such a condition depends only on the geometry of Γ_f (number of fault segments) and completely characterizes the stability. In order to perform a stability analysis let us introduce the eigenvalue problem corresponding to the static case: find φ:R × R₊ → R and β such that ∫^+∞_−∞ ∫^+∞₀ϕ²(x₁, x₂)dx₁dx₂ = 1 and

This problem has a sequence of positive eigenvalues 0<β₀ < β₁ < … with lim_n→∞β_n = +∞. The eigenvalues β_k are closely related to those of the dynamic eigenvalue problem (23)–(25). Indeed, they correspond to the intersection points of the increasing curves β→λ_k²(β) with the axis λ² = 0 (i.e. λ_k²(β_k) = 0). The first eigenvalue β₀ has a major signification in the static stability analysis:

3 The Finite Element Approach

In order to use a FEM, the finite fault zone is embedded in a bounded elastic domain Ω=]−L, L[×]0, L[. The infinite elastic half space is limited by a fictitious boundary all over which the displacement is negligible. So that we impose a null displacement along Γ_d, the part of the boundary of Ω which is not on the fault Γ_f (see Fig. 1). Let us give first the variational formulation of (23)–(25) : find Φ:Ω→R and λ² such that Φ = 0 on Γ_d and

for all functions v∈V_h (V_h is a finite element space of dimension N, composed of continuous and affine functions over each triangle) such that v = 0 on Γ_d. The domain Ω has a polygonal boundary Γ_d∪Γ_f. Therefore it is possible to cover exactly Ω with triangles. Ω is meshed using Delaunay conditions. The size of the elements decreases as we get closer to Γ_f, where a strong variation of the stress is expected. Let {e₁, e₂, …, e_N} be the canonical base of V_h and let us denote by

mass, stiffness and boundary matrices. If we put

then the eigenvalue problems (23)–(25) and (30)–(32) can be written as :

The stiffness matrix is factorized using a Cholevsky decomposition. The results are fed to a numerical solver based upon an Arnoldi-Lanczos algorithm that provides good approximations of the couples (λ_i², Φ_i). The positive eigenvalues contribute to the dominant part of the slip evolution. The largest one, λ₀², governs the essential of the slip evolution and controls the duration of the initiation phase.

4 The case of one homogeneous fault

We consider the case of one homogeneous fault of length 2a embedded in an elastic half space Ω=[−L, L]×[0, L], as shown in Fig. (1). We have chosen L large enough (L = 10a) to reduce the influence of the fictitious boundary conditions.

In order to verify the accuracy of the FEM we computed the static spectral problem. We compare the eigenvalues obtained by a FEM with 250 points on the fault with the results obtained by Dascalu et al. (2000) with a more accurate integral equation method with a finer grid of 1000 points along the fault. We find that the finite element method employed here provides good approximations of the {β_i}. The resolution of this spectral problem by two independent numerical methods validates both of them.

Fig. (2) presents the normalized eigenfunction ϕ₀(x, y) corresponding to β₀. The maximum amplitude is attained at the centre of the fault. Outside the fault (y > 0), the amplitude of ϕ₀ decays rapidly. The shape of ϕ₀ agrees with the unstable evolution of the slip velocity on the fault, as computed by a finite difference method.

4.1 Influence of the weakening parameter on the spectrum

For a finite fault of length 2a, the unstable behaviour is promoted if the weakening parameter β = aα>β₀. In the following, we consider that this condition is always satisfied. The goal is now to compute the dynamic positive eigenvalues associated with the unstable behaviour. We focus here on the dependence of the positive eigenvalues {λ_i}_{i = 1,n} on the weakening parameter β. The curves λ₀ = Λ₀(β), λ₁ = Λ₁(β) and λ₂ = Λ₂(β) are presented in Fig. 3. As expected, the number of positive eigenvalues depends on the value of β. When β<β₀, there is no positive eigenvalue. We have seen that in such a case, the fault is stable. When β₀<β<β>₁, we have only one positive eigenvalue λ₁, we have only one positive eigenvalue and so on. λ₀², the largest positive dynamic eigenvalue plays a significant role in the unstable behaviour. It is associated with the eigenfunction Φ₀(x, y) that corresponds to the fundamental mode of deformation of the system. The curve λ₀ = Λ₀(β) is monotonically increasing. This means that large values of β lead to very unstable behaviour. The asymptotic behaviour of the function Λ₀ when β→+∞ is λ₀ = β, also obtained for an infinite fault. For high weakening rate, a finite fault behaves like an infinite fault. That is, the fault finiteness is not important in the qualitative behaviour of the system. Eq. (4) leads to:

Consequently, we have:

Practically, when β is sufficiently high, we have l_c≪a. A ‘free’ initiation process is possible, like on an infinite fault. The same qualitative behaviour would be expected if α was kept constant and a→+∞. When a is large enough, the relation l_c≪a is also satisfied, and the finite fault behaves like a infinite fault.

4.2 Accuracy of the spectral approach

In the previous section, we focused on λ₀², the largest positive eigenvalue. The underlying idea is that λ₀² and Φ₀(x, y) are associated with the dominant part of the unstable evolution. The question is to know whether this assumption is correct, or not. To clear this point, we compute the slip rate evolution on a finite fault, either with a finite difference method or a finite element method. We consider a linear slip dependent friction :

where u is the slip, μ_s and μ_d the friction coefficients and L_c the critical slip. We take μ_s = 0.8, μ_d = 0.72, L_c = 0.17 m. The weakening rate α is given by eq. (2). S = ρgz is the normal stress computed at a depth of 5000 m, and G = ρc² is the rigidity modulus. The fault half length a is 500 m. We finally get: β = aα = 1.3. According to Dascalu et al. (2000), we expect only one positive eigenvalue λ₀. Fig. (4a) presents the slip rate evolution on a finite fault, computed by a finite difference method over a period of 2 s. The initial perturbation propagates on the fault. After a few reflections, another process appears. This is the development of the unstable behaviour, characterized by an exponential growth with time. Fig. (4b) presents the same slip rate evolution, computed with the finite element spectral method. The agreement between the two methods is good. The main difference is that the dominant part, characterized by (λ₀²,Φ₀(x, y)) does not take into account the propagative part of the solution. However, the dominant part gives a correct description of the slip rate distribution. We can conclude that the dominant part represents essentially the unstable behaviour of the fault, and that the assumption that we made is correct. The finite element method allows to investigate the fault behaviour when β→β₀, that is when λ₀→0. For these values, long initiation duration are expected. Fig. (5) presents the slip evolution on a finite fault with a weakening parameter β = β₀ + ε. The computation gives λ₀ = 5.10⁻⁴. The corresponding critical time (initiation duration) is T_c = 2000 s. For such a value of T_c, the finite difference method is inefficient. The slip evolution is so slow that numerical noise dominates over the initiation process. The model and the numerical method that we present here is able to produce long initiation duration (more than 30 minutes in this example). The spectral analysis that we propose is valid for any λ₀>0, even for really small values close to zero. So it is theoretically possible to produce a very broad range of initiation duration, even with this extremely simple model of a finite fault with linear friction.

5 Fault interaction

Up to now, we have considered only a simple homogeneous fault. However, it is obvious that faults are heterogeneous, composed of segments that may or may not interact. This last point is of importance in earthquake studies. The interaction of fault segments is a key issue in the prediction of the highest possible magnitude of the next event. In this section, we study the possible interaction of fault segments during the initiation process. This problem is completely different from the possible interaction of fault segments during the propagation of the rupture front.

5.1 Interaction of two fault segments

We consider here the case of a fault zone composed of two identical frictional segments (each of length 2a) separated by a rigid barrier of variable size. The goal is to show the reality of the interaction between the two segments during the eventual initiation process. We use the numerical approach already described to investigate the behaviour of the fault system. The distance between the two segments is variable. We solve the static spectral problem and compute the value of β₀ for different spacing between the segments. The results are presented in Fig. (6). As physically expected, the constant β₀ for the system decreases as the distance between the two segments decreases, which proves the interaction between the two segments. On the contrary, when the distance between the segments increases, β₀ also increases up to the constant of stability of an homogeneous fault. That is the two segments do not interact anymore. The maximum distance of interaction is of the order of 10a, where a is the half length of a segment.

5.2 The case of an heterogeneous fault

We consider a fault composed of several (11) segments, all identical, with the same weakening rate β_local, separated by rigid barriers. We aim to show that the dynamic behaviour of this complex fault system is equivalent to the dynamic behaviour of a simple homogeneous fault with a weakening rate β_equivalent. The static spectral problem associated to this complex geometry is solved. We get β₀ = 10.3. The weakening parameter of each individual segment has to be greater than β₀ to promote the unstable behaviour of the fault system. We choose β_local = 13.08. Then, the dynamic spectral problem is solved for this particular geometry.

Fig. (7) presents the first dynamic eigenfunction Φ₀^h corresponding to the heterogeneous fault system. This eigenfunction has to be compared to ϕ₀, computed for one fault segment and presented in Fig. (2). At the first glance, these two functions look very similar. However, Φ₀^h is highly heterogeneous in the close vicinity of the fault. On the other hand, the global shape (envelope) is similar to ϕ₀, which seems to indicate that the heterogeneous system behaves like an homogeneous one. The largest eigenvalue is computed to be λ₀^h ≈4.6 and is associated with Φ₀^h. With the help of the curve λ₀ =Λ₀(β) (defined for one homogeneous fault segment) presented in Fig. (3), it is possible to derive the value of the corresponding weakening parameter:

A homogeneous fault with such a weakening rate will lead to λ₀ = λ₀^h. We argue that this particular value of β is the equivalent weakening rate. That is:

This last hypothesis will be discussed in next section.

6 The initiation pattern: a possible signature of instability

6.1 Case of one homogeneous fault

Characterizing the unstable behaviour of a fault will be a step in earthquake prediction. Through eq. (1) we have seen the simple dependency of the dominant part of the unstable slip evolution with respect to the y coordinate in the case of an infinite fault. This can be presented as follows

i.e. during the initiation phase the linear weakening condition on the fault is transmitted everywhere in the interior of the elastic medium. This property can be a candidate for the signature of the instable evolution of the finite fault model. The question is how to find a similar property in case of a finite fault. Let us suppose that we deal with a slow initiation (i.e. N = 1, λ₁² < 0 < λ₀²≪1) or we can neglect the contribution of all the other eigenfunctions with respect to the first one (i.e. 0<λ²_N−1<…<λ₁²≪λ₀²). In this case

and the associated ratio given by eq. (38) can be deduced as:

This ratio represents the information about the weakening rate of the fault, when it is defined. Actually, B_β is defined everywhere on Ω−Γ_d. On Γ_c (x₂ = 0) we have the relation B_β = β. Fig. (8) presents the function B_β0(x₁, x₂) computed in the static case. We have focused on the behaviour of B_β0 close to the fault Γ_c. The minimum value of the function is β₀. It is reached on the fault and on a small domain close to the fault. Note the typical shape of the function, organized in a band of the size of the fault length. Fig. (9) presents some computations for different weakening rates β for 1.2–5. The general shape of B_β is the same for all β and then defines an initiation pattern that qualitatively characterizes the unstable behaviour of the fault. The most interesting point is the existence of a domain of B over which B_β(x, y) = β. This signifies that the weakening parameter can be measured in the surrounding elastic body. We can define now the domain of confidence as the part of the elastic body where the previous property is satisfied:

The fault surface Γ_f is always included in D_c. The extension of the domain of confidence is an increasing function of the weakening rate β. From the relation between λ₀ and β (see Fig. 3) and the measure of β over the domain of confidence, it is possible to get the largest positive eigenvalue. Since the initiation duration is inversely proportional to λ₀, we are able to prescribe a critical time, that is we can predict the time of occurrence of the future event.

6.2 Case of an heterogeneous fault

We consider the case of a fault zone composed of eleven identical frictional segments (described in Section 5.2). We now compute the function B(x₁, x₂) for this complex fault system. The result is presented in Fig. (10). As already seen, close to each individual fault segment, an initiation pattern develops in the elastic medium, associated with a local domain of confidence over which β_local can be measured. The striking feature of Fig. (10) is the existence of a wide domain over which B(x₁, x₂) is nearly constant, independent from the individual fault segments but closely related to the whole fault system. All over this domain, we have the following relation:

That is, all over this wide domain, it is possible to measure the collective behaviour of all the fault segments, similar to the behaviour of an homogeneous fault with β = β_equivalent. Some strong results arise from these computations:

• (i) An initiation pattern develops inside the elastic medium and is characteristic of the unstable behaviour of the fault system.

• (ii) The collective behaviour of all the segments is transmitted in the medium and can be measured over a wide domain.

• (iii) There is a spectral equivalence between a complex heterogeneous fault (β_local, λ₀^h, Φ₀^h) and an homogeneous fault (β_equivalent = Λ₀⁻¹(λ₀^h), λ₀ = λ₀^h, Φ₀).

These results are in perfect compliance with those of Campillo et al. (2001), who used finite difference time domain computations to demonstrate the accuracy of an effective friction law deduced from the properties of the strain field inside the elastic body.

7 Discussion

The unstable behaviour of a fault is successfully modelled by the dominant part of the slip rate evolution. The model that we develop in this paper is rather simple: a finite homogeneous fault under a linear slip dependent friction law, characterized by the weakening rate α and its half length a through the parameter β = aα. This model produces long initiation duration, up to at least 30 minutes, as presented in Fig. (5). Many other processes may be involved in the initiation process over such a long period of time. We do not argue that these physical mechanisms, such as fluids effects, are not important, we just show that a simple fault model can produce a very broad range of initiation durations with no need for other mechanisms.

We have shown that the analysis we performed on a finite homogeneous fault remains valid for a heterogeneous fault. The collective behaviour of interacting fault segments is very similar to the homogeneous finite fault behaviour. The equivalent weakening rate is given by the relation λ₀ = Λ₀(β). The equivalence is complete in terms of initiation process, that is λ₀ is the same for the two faults. Therefore, the critical time (the initiation duration) is the same in both cases. The equivalent homogeneous fault is a sort of homogenization of the complex fault. This technique may be used for different scales of heterogeneity to provide a correct description of the global behaviour of complex, highly heterogeneous fault systems. However, the homogenization implicates the loss of information on how the rupture process will develop in detail. The rupture complexity observed both in laboratory experiments or in strong motion inversions will be lost by these homogenization techniques.

The transmission of information about the weakening rate the elastic medium offers the possibility to constrain the time of occurrence of the next event, even without knowing anything on the fault. The existence of an initiation pattern associated with the unstable behaviour of the fault and the measure of the function B(x, y) is sufficient to give a good approximation of λ₀, and therefore to prescribe the critical time. However interesting these results may be, they have to be confirmed by other methods, and moreover they have to be extended to the 3D case. Many other parameters have to be taken in account, such as the fault geometry or the orientation of the fault plane. This study appears like a first attempt to recognize and characterize a possible signature of the unstable behaviour of a fault, that is an initiation pattern.

8 Conclusions

We presented a numerical method based on the finite elements techniques to study the problem of the initiation process on a finite fault. This method allows accurate computations for large eigenvalues which is not the case for the numerical approach based on the integral equations proposed by Dascalu et al. (2000). In case of an unstable fault, we show that the catastrophic evolution of the slip rate is successfully modelled by the dominant part, formed by the contributions of the positive eigenvalues. We give the relation λ₀ = Λ₀(β) between the weakening rate β and the largest positive eigenvalue λ₀ that largely controls the unstable behaviour. We show that some information about the weakening rate β is transmitted in the elastic medium over a domain of confidence. We present the initiation pattern associated with the unstable behaviour. This analysis remains valid for more complex fault systems. We show the complete equivalence, in terms of initiation process, between a heterogeneous fault characterized by (β_local, λ₀^h, Φ₀^h) and a homogeneous fault characterized by (β_equivalent, λ₀^h, Φ₀) where β_equivalent =Λ₀⁻¹(λ₀^h). These preliminary results have to be confirmed by other methods, and moreover, they have to be extended to the 3D case.

Acknowledgments

This work is financially supported by the CNRS and ANDRA through the GdR FORPRO (Research action number 99.VI) and corresponds to the GDR FORPRO contribution number 2000/30 A. This paper is dedicated to Marine, William and their parents, Geoffrey and Karine.

References

Author notes