Abstract. We present an efficient and rigorous numerical procedure for calculating the elastodynamic response of a fault subjected to slow tectonic loading processes of long duration within which them are episodes of rapid earthquake failure. This is done for a general class of rate-and state-dependent friction laws with positive direct velocity effect. The algorithm allows us to treat accurately, within a single computational procedure, loading intervals of thousands of years and to calculate, for each earthquake episode, initially aseismic accelerating slip prior to dynamic rupture, the rupture propagation itself, rapid post seismic deformation which follows, and also ongoing creep slippage throughout the loading period in velocity-strengthening fault regions. The methodology is presented using the two-dimensional (2-D) antiplane spectral formulation and can be readily extended to the 2-D in-plane and 3-D spectral formulations and, with certain modifications, to the space-time boundary integral formulations as well as to their discretized development using finite difference or finite element methods. The methodology can be used to address a number of important issues, such as fault operation under low overall stress, interaction of dynamic rupture propagation with pore pressure development, patterns of rupture propagation in events nucleated naturally as a part of a sequence, the earthquake nucleation process, earthquake sequences on faults with heterogeneous frictional properties and/or normal stress, and others. The procedure is illustrated for a 2-D crustal strike-slip fault model with depth-variable properties. For lower values of the state-evolution distance of the friction law, small events appear. The nucleation phases of the small and large events are very similar, suggesting that the size of an event is determined by the conditions on the fault segments the event is propagating into rather than by the nucleation process itself. We demonstrate the importance of incorporating slow tectonic loading with elastodynamics by evaluating two simplified approaches, one with the slow tectonic loading but no wave effects and the other with all dynamic effects included but much higher loading rate.

Guided by seismic observations of short-duration radiated pulses in earthquake ruptures, Heaton (1990) has postulated a mechanism for the frictional sliding of two identical elastic solids that consists in the subsonic propagation of a self-healing slip velocity pulse of finite duration along the interface. The same type of pulse may be conjectured for inhomogeneous slip along sufficiently large, and compliant. technological surfaces. We analyze such pulses, first as steady traveling waves which move at constant speed, and without alteration of shape, on the interface between joined elastic half-spaces, and later as transient disturbances along such an interface, arising as slip rupture propagates spontaneously from an over-stressed nucleation site. The study is conducted in the framework of antiplane elastodynamics ; normal stress is uniform and alteration of it is not considered. We show that not all constitutive models allow for steady traveling wave pulses: the static friction threshold subsequent to the relocking of the fault must increase with time. That is, such solutions do not exist for pure velocity-dependent constitutive models, in which the stress-resisting slip on the ruptured surface is a continuously decreasing function of the instantaneous sliding rate (but not of its previous history or of other measures of the evolving state of the surface). Further, even for constitutive models that include both the rate-and state-dependence of friction. such as the laboratory-based constitutive models for friction as developed by Dieterich (1979Dieterich ( , 1981 and Ruina (1983), steady pulse solutions do not exist for versions, like one discussed by Ruina (1983), which do not allow (rapid) restrengthening in truly stationary contact. For a particular class of rate-and statedependent laws which includes such restrengthening, we establish parameter ranges for which steady pulse solutions exist, and use a numerical method stabilized by a Tikhonov-style regularization to construct the solutions. The numerical method used for the transient analysis adopts Fourier series representations for the spatial dependence of stress and slip along the interface, with the (time-dependent) coefficients in those Fourier series being related to one another in a way which obtains from exact solution to the equations of elastodynamics.This allows an efficient numerical method, based on use of the Fast Fourier Transform in each time step, with the frictional constitutive law enforced at the FFT sample points along the interface. Solutions based on a law that includes restrengthening in stationary contact show that spontaneous rupture propagation will occur either in the self-healing slip pulse mode (but not generally as a steady pulse) or in the classical enlarging-crack mode, depending on the values of parameters which enter the constitutive law. This analysis suggests that the strictly steady, traveling wave pulse solutions may either be unstable or have a limited basin of attraction.

We use the finite element method to analyze stress variations in and near a strongly coupled subduction zone during an earthquake cycle. Deformation is assumed to be uniform along strike (plane strain on a cross section normal to the trench axis), and periodic earthquake slip is imposed consistent with the long-term rate of plate convergence and degree of coupling. Simulations of stress and displacement rate fields represent periodic fluctuations in time superimposed on an average field. The oceanic plate, descending slab, and continental lithosphere are assumed here to respond elastically to these fluctuations, and the remaining mantle under and between plates is assumed to respond as Maxwell viscoelastic. In the first part of the analysis we find that computed stress fluctuations in space and time are generally consistent with observed earthquake mechanism variations with time since a gmat thrust event. In particular, trench-normal extensional earthquakes tend to occur early in the earthquake cycle toward the outer rise but occur more abundantly late in the cycle in the subducting slab downdip of the main thrust zone. Compressional earthquakes, when they occur at all, have the opposite pattern. Our results suggest also that the actual timing of extensional outer rise events is controlled by the theology of the shallow aseismic portion of the thrust interface. The second part of the analysis shows the effects of mantle relaxation on the rate of ground surface deformation during the earthquake cycle. Models without relaxation predict a strong overall compressional strain rate in the continental plate above the main thrust zone, with the strain rate constant between mainshocks. However with significant relaxation present, a localized region of unusually low compressional, or even slightly extensional, strain rate develops along the surface of the continental plate above and somewhat inland from the downdip edge of the locked main thrust zone. The low strain rate starts in the middle or late part of the cycle, depending on position. This result suggests that the negligible or small contraction measured on the Shumagin Islands, Alaska, during 1980 to 1991, may not invalidate an interpretation of that region as being a moderately coupled subduction zone. In contrast, mantle relaxation causes only modest temporal nonuniformity of uplift rates in the overriding plate and of extensional stress rates in the subducting plate, even when the Maxwell time is an order of magnitude less than the recurrence interval.

Studies of the mechanics of subduction as inferred from earthquake cycle observations suggest that the distribution and style of seismicity in the seafloor, between the trench and the outer rise, and in the slab at intermediate depth, can in some cases serve to identify asperity locations along the thrust interface [Dmowska and Lovison, 1992]. Such asperities, identified from seismic wavefield modeling, are the zones of highest seismic moment release in large underthrusting events. To the extent that asperity locations are relatively stationary from one event to the next, their locations provide the zones of highest expected moment release in future large earthquakes, and rupture often nucleates at the border of an asperity. The region of the thrust interface outside such asperities is, apparently, less well coupled and releases moment throughout the great earthquake cycle in some combination of aseismic creep and moderate seismicity. Thus it is reasonable that stress and deformation rates associated with the earthquake cycle should be most pronounced near asperities, and that this should have seismic and geodetic consequences. Three‐dimensional finite element modeling is used here to understand such stress and deformation patterns and their variation in time, in relation to heterogeneity of coupling along thrust interfaces. The stress field helps to explain the observed clustering of seafloor seismicity along the strike of the convergent margin. In cases of convergence at approximately normal incidence, like for the region of the Valparaiso, Chile, 1985 thrust event, the modeling is consistent with the observation that areas of large earthquakes in the seafloor toward the outer rise and in the slab tend to lie within corridors through thrust zone asperities, running perpendicular to the line of the trench. We seek to learn if such model stress fields are consistent with observations, for the strongly oblique subduction margin of the Rat Islands, western Aleutians, 1965 event, that active areas of the outer rise and slab at intermediate depth are offset along strike from asperity locations. Modeling results here for the stress in the seafloor raise the possibility that to explain this offset, the asperity zones along the thrust interface may have to be strung out along the direction of oblique slip, perhaps reflecting the contact path of subducting seamounts or geometric irregularities along the interface. Shear stress patterns created in the upper plate, when there is oblique subduction, suggest that favorable areas for back‐arc strike slip activity following underthrusting, as in the Adak Island, central Aleutians, region of the 1986 Andreanof Island earthquake [Ekström and Engdahl, 1989], will also be shifted along strike from asperity locations. Our analyses show how deformation patterns on the earth's surface above asperities differ from patterns above nonasperities, and hence provide tools to identify inhomogeneous coupling from geodetic observations. We discuss possible bathymetric, topographic, and struct...

Deformation associated with the earthquake cycle in the Shumagin Islands segment of the Alaska-Aleutian subduction zone is analyzed with the use of a two-dimensional finite element model. The model consists of an oceanic plate dipping under an upper plate, both of which respond elastically to stress fluctuations in the earthquake cycle, and these are underlain by asthenospheric manfie and manfie wedge regions which respond viscoelastically. It is tailored to the geometry of the Shumagin Islands region, by using seismicity to define the position of the interplate interface and (partially) coupled region along it. The model is preconditioned by forcing this interface to undergo periodically repeated slips up to (and including) the time of the May 31, 1917, event (Ms= 7.4) in that region, with each chosen to be consistent with the moment and estimated rupture area of that event. We investigated the dependence of model results for geodetic signals on the strength of seismic coupling between the plates and viscoelastic relaxation of deviatoric stresses in the mantle, including in the mantle wedge close to the plate junction and along the aseismic downdip continuation of the thrust interface. In models with significant relaxation in the wedge or downdip thrust zone, results show that as the intraseismic stage matures, there is a region of diminished compressional strain rates, and even of locally extensional rates, on the Earth's surface above the downdip end of the seismically coupled zone. Based on the seismic estimates of the location of the coupled zone, this region is in the area of the Shumagin Islands. We find that if approximately 20% of the convergence takes place seismically (compatible with the previous seismic history), and if an extensive region of relaxed deviatoric stress is assumed to be present in the wedge and/or along the downdip interface, then deformations predicted by the model can be made consistent with the measured strain data from the Shumagin Islands geodetic network [Lisowski et al., 1988; Larson and Lisowski, 1994], as well as uplift and tilt data [Savage and Plafker, 1991; Beavan, 1992]. Our model simulations here suggest that the Shumagin segment is capable of large earthquakes. The hypothesis of totally aseismic subruction is not similarly consistent with all geodetic constraints. of a locked gap [e.g., Sykes and Jaumd, 1990; Dmowska and Paper number 95JB03461 0148-0227/96/95JB-03461 $05.00 Lovison-Golob, 1991; Jaumd and Estabrook, 1992; Bufe et al., 1994]. The most recent event, continuing the trend of accelerated moment release that had already been reported, is the M s = 6.9, M o = 0.2-0.3 x 1020 N m event of May 13, 1993 [Lu et al., 1994; Tanioka et al., 1994]. Also, since 1977, there has been a cessation of the previously abundant, generally tensional, outer-rise seismicity [Dmowska and Lovison-Golob, 1991]. We have interpreted this as a sign that significant thrust zone coupling is decreasing the extensional stresses in the outer rise as the seismic cycle matures [Dmowska et al., 1...