friction law will be discussed after we set up the elastodynamical problem.

We have assumed antiplane slip, so that the only component of displacement that we consider is *u=u*_{y} and the relevant traction component across the fault is *T=−σ*_{yz}*.* Thus, *u* satisfies the wave equation

**[1]**

where *β* is the shear wave speed. By using the symmetry about the fault plane *z*=0, we find the appropriate boundary conditions for the solution of the problem in the upper half-space (*z*>0):

*u*(*x,* 0^{+}, *t*)=0 on -∞ *<x<−ℓ*_{1}*(t)* and *ℓ*_{2}*(t)*<*x*<∞ **[2]**

outside the slipping part of the fault delimited by *ℓ*_{1}*(t)* and *ℓ*_{2}*(t),*

*T*(*x,* 0^{+}, *t*)*=T*_{e}*(t)+T(D, V, θ)* on −*ℓ*_{1}*(t) <x<ℓ*

inside the fault. *T(D, V, θ)* is the friction law. This problem is then a typical crack problem with mixed boundary conditions, which consists in finding the slip *D*(*x,* 0, *t*)*=*2*u*(*x,* 0^{+}, *t*) on the fault.

In spite of the apparent simplicity of this crack problem, in most cases it is intractable by analytical methods because of the mixed boundary conditions (Eqs. **2** and **3**). As discussed (25), finite differences are too inaccurate because of high frequency dispersion, so that we adopted a boundary integral equation method that we believe is the most appropriate for solving this kind of crack problem (26, 27). The appropriate boundary integral equation (10) is

**[4]**

where *S* is the space-time convolution integral

**[5]**

with . The convolution over space has a Cauchy-like singularity at *ξ*=*x.* The convolution over time, on the other hand, is regular because *K(x, t)* → 0 when *x* → 0. The details of the numerical method used to solve this problem and its implementation in a CM5 parallel computer have been discussed (10). Thus it is sufficient to say that Eq. **4** is a nonlinear algebraic equation for the simultaneous solution of *D, V, T,* and *θ*.

Let us remark that the fault is assumed to be clamped at the ends by unbreakable barriers. This boundary condition is different from the periodic conditions used by Carlson and Langer (14). To apply periodic boundary conditions, we would have to change the Green functions used to derive Eq. **4.**

The nature of the solution of Eq. **4** is completely determined by the friction law *T(D, V, θ)*. Extensive work on rock friction has been discussed by Dieterich (28) and does not need to be repeated here. Experimental evidence (20, 21, 29) shows that friction laws at low slip rates should at least include three elements: *(i)* direct stress change for rapid increases in slip velocity, *(ii)* an intrinsic time constant for the response to abrupt changes in velocity, and *(iii)* velocity weakening at steady-state slip.

At present we are unable to introduce such laws in our numerical method because they were proposed for very low values of slip velocity, so that the solution of the boundary integral equation problem would have to take into account slip rates varying over more than six orders of magnitude (i.e., from *μ*m/s for fault creep up to m/s for seismic slip). Since the time step for the solution of Eq. **4** is controlled by the faster time scales present in it, the slow evolution of the fault in the interseismic period would have to be computed at the time steps appropriate for the dynamic regime, which is several orders of magnitude less than those that are actually needed. Equations of this type are called stiff and require special techniques for their solution as shown by Tse and Rice (30) for a spring-loaded massive slider. A possible method to increase the time steps during the interseismic period would be to modify Eq. **4** at low speeds by using its quasi-static approximation as proposed by Perrin *et al.* (31), but we have not implemented this option yet.

Because of the difficulties discussed in the previous paragraph, we have used a simplified friction law. The most sweeping assumption is that the fault is perfectly locked during the inter-seismic period so that successive events can be treated independently. In the Dieterich-Ruina rate- and state-dependent laws, slip occurs at all times no matter how small the total stress. In the period when the fault is locked, we assume that there is no slip on the fault, so that stress increases steadily during the continuous load *T*_{e}*(t)*. When the total stress on a certain point of the fault reaches a threshold *T*_{u}*,* slip begins and we immediately switch to the dynamic equation (Eq. **4**).

The particular friction law assumed in our computations is

*T(D, V, θ)*=*T*_{u}(1*−D/U*_{0}) for *D*<*U*_{0}

*T(D, V, θ)*=0 for *D>U*_{0} and *θ>θ*_{0}**[6]**

*T(D, V, θ)*=*T*_{sp}(1*−θ*/*θ*_{0}) for *D>U*_{0} and *θ<θ*_{0},

to which we add the following evolution equation for *θ*

**[7]**

where *D* is slip measured from the beginning of the current slip episode, *θ* is a state variable, and *V* is instantaneous slip velocity. At steady state, *θ* is equal to the slip rate, so that *θ* differs from slip rate *V* only during fast rate changes. The first-order evolution equation for *θ* introduces a relaxation time scale *τ*_{0} and a corresponding relaxation length *D*_{c}*=βτ*_{0}. These are the small scale parameters needed to regularize the fault model at small scales and avoid the intrinsic discreteness of our previous models (10). The dependence of friction on slip and steady-state slip rate *(θ)* is better appreciated in Fig. 2. As slip begins, there is a slip-weakening process that simulates