A Study of Interactions Between Thrust and Strike-slip Faults

A 3-D finite difference method is applied in this study to investigate a spontaneous rupture within a fault system which includes a primary thrust fault and two strike-slip sub-faults. With the occurrence of a rupture on a fault, the rupture condition follows Coulomb’s friction law wherein the stress-slip obeys the slip-weakening fracture criteria. To overcome the geometrical complexity of such a system, a finite difference method is encoded in two different coordinate systems; then, the calculated displacements are connected between the two systems using a 2-D interpolation technique. The rupture is initiated at the center of the main fault under the compression of regional tectonic stresses and then propagates to the boundaries whereby the main fault rupture triggers two strike-slip sub-faults. Simulation results suggest that the triggering of two sub-faults is attributed to two primary factors, regional tectonic stresses and the relative distances between the two sub-faults and the main fault.


INTRODUCTION
The primary fault in a highly complex fault system is usually accompanied by several branch faults.The 1992 Landers earthquake (M w 7.3), 1998 Antarctic earthquake (M w 8.1) and 1997 Kagoshima earthquake (M w 6.0) are all good examples (Kase and Kuge 2001).It is necessary to understand the interactions between these faults when an earthquake occurs.For this purpose, the rupture processes of the main fault and those of the branch faults caused by stress changes induced from the failure of the main fault must be investigated.Meanwhile, the dynamic behavior of a fault system within a short time period, from a few seconds to just tens of seconds, should be inferred from observations.The rupture processes include initiation, propagation, and arrest of ruptures (Aki 1979;Lindh and Boore 1981;Scholz 1990).Entire rupture processes are controlled primarily by physical and geological conditions of the fault system.
Before the 1990s, fault braches were not included in the models used for the studies of dynamic rupture processes (Das and Aki 1977;Miyatake 1980;Day 1982;Mi-kumo et al. 1987).But current studies had considered the occurrence of a mainshock and additional stresses which are generated in the source area and its surroundings and, perhaps, trigger subsequent events within branch faults (King and Cocco 2001;Wang and Chen 2001;Wang et al. 2003Wang et al. , 2010;;Lin and Stein 2004;Ma et al. 2005;Chan and Stein 2009).Recently, this type of fault interaction, especially for a strike-slip fault system has been widely studied by earth scientists.One of the notable examples is the San Andreas fault in California (Harris and Day 1993;Yamashita and Umeda 1994;Kase andKuge 1998, 2001).Deffontaines et al. (1994Deffontaines et al. ( , 1997) ) first pointed out possible interactions between a thrust front and transfer faults in Taiwan.
The stress changes in a source area and its surrounding areas due to a failure of a thrust fault has garnered much attention by seismologists.For example, Ma et al. (1999) studied the spatial stress re-distribution caused by the devastating 1999 M w 7.6 Chi-Chi earthquake which ruptured along the thrust-type Chelungpu fault.From the focal mechanism of the mainshock and aftershocks of the Chi-Chi earthquake, Mozziconacci et al. (2009) investigated the interaction between the Chelungpu fault and the Sanyi-Puli transfer fault zone (a strike-slip fault zone).Mozziconacci et al. (2013) investigated the interaction between the thrust main fault of the 2006 M w 6.1 Taitung earthquake and its strike-slip sub-faults.The geologic map of active faults in Taiwan (Lin et al. 2000) shows several strike-slip faults (from south to north) which are indeed located near respective main thrust faults.For example, in both the northern and southern ends of the Chelungpu fault there are two sub-faults (Lee et al. 2003(Lee et al. , 2011)).The two sub-faults are not purely strike-slip faulting because both exhibit a dip-faulting component.
To simplify the problem, we investigate the interaction between a thrust and two strike-slip faults.In the model, the thrust fault and two strike-slip faults are perpendicular to one another, and the fault geometry is similar to the Chelungpu fault system.This type of study is particularly important for seismically active areas which have a higher population density, are highly industrialized, and have a history of destructive earthquakes.Taiwan is a typical example (Wang 1998(Wang , 2008)).The 1906 M 7.1 MeiSan earthquakes (cf.Wang 1998), the 1935 M 7 HsinChu-Taichung earthquake (Miyamura 1985;Wang 1998), and the 1999 M w 7.6 Chi-Chi earthquake (Ma et al. 1999) all occurred in highlypopulated areas of Taiwan.
As mentioned earlier, the study of rupture processes in a complex fault system is not only necessary for academic purposes but also for seismic hazard estimates.In this study, a 3-D finite difference method is applied to simulate fault interaction and rupture processes in a model consisting of a main thrust fault linked with two strike-slip faults.The numerical method developed by Wang (2005) will be applied in this study and described in detail below.

METHODOLOGY
The present model consists of a main thrust fault and two strike-slip faults.Since we are unable to fully elucidate the free surface effects on rupture propagation, the three faults are located 4 km below the free surface in a homogeneous isotropic medium as shown in Fig. 1a.The main thrust fault (i.e., the central dipping plane in Fig. 1a and named the M-fault hereafter) is defined in the fault coordinates represented by the X Y Z l l l coordinates.For simplification, this system is called the F-coordinate.If the thrust M-fault has a dip angle of 45°, the grid points of the top edges of the two strike-slip type sub-faults are collocated with the grid points of the F-coordinate system and only one coordinate system can be used to describe both M-and strike-slip faults.Since the M-fault has a dip angle of 30°, the grid points of the top edges of the two strike-slip-type sub-faults cannot be located exactly at the F-coordinate grid points (see Fig. 2, along the Xl-axis).We need to define a new coordinate system.The two vertical strike-slip sub-faults (i.e., the slashed planes in Fig. 1a) that are perpendicular to the M-fault in the XYZ coordinates are the original coordinates called the O-coordinate.The northern and southern sub-faults (named the N-fault and the S-fault, respectively) are not directly connected to the M-fault, and each separates from the M-fault with 1 km.In the following simulations, the dip angle, δ, of the M-fault is set at 30° and the strike is 0° (see Fig. 2) which is almost the average dip angle of numerous major thrusts in Taiwan such as the Chelungpu fault.The fault length and width of the M-fault are 8 and 4 km, respectively in the strike and dip directions.From the relationships of earthquake magnitude versus fault length inferred by Well and Coppersmith (1994), the magnitude of the simulated event is 5.7 because the fault length of the main thrust fault is 8 km.
The two sub-faults have the same length and width, which are 3.46 and 2 km, respectively.The depths of the top edges of the three faults are the same.The initial external regional stresses as shown in Fig. 1b are denoted by XX v , YY v , and ZZ v .Meanwhile, XX v and YY v are the maximum compressive stress and the second largest stress.The angle between the strike of the M-fault and the maximum compressive stress is denoted by φ.Note that θ and φ are a pair of complementary angles.While ZZ v is the vertical stress normal to the M-fault.
The second-order 3-D finite difference method is applied to simulate ruptures on the M-fault and the transfer of ruptures to the N-and S-faults (Kase andKuge 1998, 2001).In order to ensure that the calculations are made with highresolution (Olsen et al. 1997;Madariaga et al. 1998) and the 30°-dipping M-fault is completely and exactly located at the grid points, the grid lengths of the M-fault are set to be dxl = 0.2 km, dyl = 0.2 km, and dzl = 0.1 km in the F-coordinate.Because the top edges of the two sub-faults must be as high as the main fault, the two sub-faults cannot be located exactly at the grid positions in the F-coordinate (see Fig. 2, along the Xl-axis), we must let the grid lengths of the sub-faults be dx = 0.173 km, dy = 0.2 km, and dz = 0.1 km in the O-coordinate.The physical quantities on the grid points of the M-fault and sub-faults are connected through a 2-D interpolation technique.The numbers of grid points are 41 × 21 (along the Yl-and Xl-axes, respectively) in the F-coordinate for the M-fault and 21 × 21 (along the X-and Z-axes, respectively) in the O-coordinate for the Nand S-faults.
When the rupture is initiated at the center of the Mfault (marked by a star in Figs.1a and b) at time t = 0, a decrease in the shear stress, v , from the maximum static frictional stress (i.e., the stress drop) results in a displacement wherein the rupture then propagates outward when t > 0. This would induce stress changes in other non-ruptured areas (including sub-faults).When the shear stress exceeds the maximum static frictional stress in the sub-faults, additional ruptures could also occur.After a rupture occurs, the variation in shear stresses follows slip-weakening fracture criteria (Andrews 1976;Day 1982), as illustrated in Fig. 3, where 0 v , s v , and d v are, respectively, the initial shear

BOUNDARY CONDITIONS
In a homogeneous and isotropic medium, the wave equation is: where u is displacement, t is density, and m and n are Lamé constants.Except at the location of the initial rupture on the fault plane; the displacements and velocities are all zero: The initial external stresses are represented by: Where XX v and YY v are the stresses along the directions of two horizontal compressive stresses which are in an inward direction as shown in Fig. 1b and ZZ v is in the vertical direction.
To calculate the stresses on the M-fault, the initial stresses must be rotated from the O-coordinate onto the Fcoordinate using the following two steps.First, the initial regional stresses are rotated φ degrees counter-clockwise to the Z-axis (Fig. 1b).Thus, Where A z z ^h is the rotational matrix with a symmetric axis in the vertical direction.Second, the initial regional stresses are rotated d degrees counter-clockwise to the Y-axis.Thus, where A y d ^h is the rotational matrix with a symmetric axis in the direction of the strike.The normal stress on the Mfault is now: and the shear stress is: where d n is the dynamic frictional coefficient.For the noninitial crack grids, the direction, r i of the rupture is identified from the shear stress from zx vm and zy vm in Eq. ( 7b) and Similarly, for the N-and S-faults the applied regional stresses are rotated by φ degrees counter-clockwise onto the O-coordinate, as shown in Eq. ( 6).For these two sub-faults, the normal stress is: and the shear stress is: The boundary conditions for the M-fault in the F-coordinate are: When the fault breaks, the boundary conditions are Eqs.(9c) and (10c), and the shear stress obeys slip-weakening fracture criteria (see Fig. 3); thus, ( ) is the slip in the direction of the rupture; u i + and u i -are the displacements on the hanging wall and the footwall of the fault, respectively; s n is the static frictional coefficient.The boundary conditions for the N-and S-faults in the original coordinates are: likewise the stresses are continuous: Similarly, when a fault rupture occurs, the boundary conditions are Eqs.(12b) and (13b); the shear stress obeys the slip-weakening fracture criteria.Thus, we have: As mentioned above, the fault edges are located 4 km below the free surface.Nevertheless, the free surface boundary condition proposed by Ilan and Loewenthal (1976) is used.The absorption boundary method is used for 3-D finite differences (Higdon 1991).The wave equations are replaced by central finite difference equations.The boundary conditions on the fault are derived following the procedures suggested by Horikawa (1996) (see the Appendix for an example of the calculations of z z v l l and u xl on the M-fault).A detailed description of the derivation of stresses and slip can be found in Wang (2005).
The rupture starts from the M-fault center (marked with a star in Fig. 1a) where the finite difference is performed on the F-coordinate (i.e., the Xl, Yl, and Zl axes) per the black circles connected with the dashed lines in the vertical profile of Fig. 2. The grid points on the O-coordinate are denoted by white circles, and the difference values on the white circles are retrieved from those at the black grid points in two hanging blocs and two foot blocs of the M-fault by means of a 2-D interpolation (Press et al. 1992).For example, the values at black grid points , , , and (in the F-coordinate system) in the hanging bloc are used to interpolate the values at white grid point A (in the O-coordinate system).Similarly, the value at grid point B in the foot bloc is interpolated from those at , , , and .Note that the values at two pairs of grid points, i.e., and as well as and , on the M-fault plane are the same before the rupture, yet become different afterwards because of the slip discontinuity on the ruptured plane.A similar interpolation is performed for all of the white grid points that are right above and right below the fault plane.At this step, the slips at the grid points are in the F-coordinate system and must be rotated into the O-coordinate system; they must be rotated onto the O-coordinate for the computations of the 3-D finite differences in these coordinates for the remaining grid points.

PARAMETERS AND STABILITY CONDITIONS
The values and orientations of initial external stresses given in Table 1 were taken from Kase andKuge (1998, 2001) and were tested through numerical simulations by Wang (2005).The stresses can trigger the strike-slip subfaults due to the break of the thrust main fault.From Coulomb's friction law (Harris and Day 1993), the shear stress and the maximum static frictional stress in spatial and temporal distribution are described by the following equation: is negative, the segment can be triggered to rupture.From the initial external stresses given in Table 1, the relation between φ and s D on the two strikeslip sub-fault planes at t = 0 s is shown in Fig. 4. Obviously, the minimum of s D is at φ = 30°, thus indicating that the sub-faults can be most easily triggered at φ = 30°.Therefore, the simulations and numerical tests of parameters in use are made only at φ ± 1° for each φ.
From the tests, we find that the two strike-slip subfaults (i.e., the N-and S-faults) can be triggered to rupture only when 18° < φ < 42° and the range of φ depends on the fault system.The initial stresses on the fault system are XX v , YY v and ZZ v .We must transfer them from the tectonic system to the fault plane.This gives the angle ( ) m , which is the maximum shear stress direction.According to the tectonics, the maximum stress direction will be a raked direction of the fractured fault.In this study, the rake angle of the main fault is set to be 60°, so that the φ angle is 26.57°.The sampling time interval, dt, is required to meet the computational stability of the finite differences (Kelly et al. 1976), i.e., dt ds where ds is the grid size; and a and b are, respectively, the P-and S-wave velocities.Thus, for ds ≥ 0.05 km dt = 0.01 s can meet the stability condition of Eq. ( 16).
Although the rupture velocity cannot be determined from this model, several authors (e.g., Andrews 1976;Das and Aki 1977;Day 1982;Harris and Day 1993) stated that the rupture velocity can be identified from the difference between the maximum static frictional stress ( s v ) and the initial stress ( 0 v ) as well as from the difference between dynamic frictional stress ( d v ) and initial stress ( 0 v ).Andrews (1976) and Das and Aki (1977) defined the parameter S by the four stresses to represent the differences on rupture velocities, i.e., ( ) ( ) When S ≥ 1.63, the rupture velocity is less than the Rayleighwave velocity, and it increases with rupture length and eventually approaches the latter.By contrast, when S < 1.63, the rupture velocity starts from the sub-Rayleigh velocity and increases with rupture length and can in fact be greater than the S-wave velocity or even approach the P-wave velocity (Aki and Richards 1980;Wang 1996Wang , 2012)).The sub-Rayleigh velocity is the velocity of rupture waves which is slower than the velocity of Rayleigh waves.In this study, for the ruptures on the thrust fault and the triggering on two sub-faults, the set of the initial external stresses is shown in Table 1.Obviously, the value of S is smaller than 1.63 on the M-fault, and thus the rupture velocity can approach the P-wave velocity.

SIMULATION RESULTS
Using the parameters listed in Table 1 along with φ = 26.57°(or θ = 63.43°), which gives a rake of 60° in the finite difference computations, and the initial rupture occurs at the center of the M-fault (designated as a star in Fig. 1).The equal-time contours of ruptures on the three faults are illustrated in Fig. 5, in which the view is set to be in the direction perpendicular to the fault plane and from the hanging wall.In Fig. 5, the numbers represent the times of ruptures in seconds.When t = 0 ~ 1 s (Figs.5a -c), the M-fault starts to rupture from its center (Fig. 5a) and when the rake is 60°, the rupture velocity is faster along the dip than along the strike.In this time period, the ruptures propagate outward, while the S-and N-faults are still locked.Overall, it takes about 1.71 s to make the M-fault completely rupture.Figure 5d shows that the S-and N-faults start to rupture at t ~1.37 s.The S-fault (see the left panel) and N-fault (see the right panel) are both triggered at a depth of 5 km, but the ruptures on the two sub-faults propagate along the opposite directions.The two sub-faults continue to rupture for about 2.23 s (Fig. 5f).It is noted that the asymmetric rupture patterns on the S-and N-faults are caused by different ruptures from the M-fault: downward and upward ruptures on the S-and N-faults, respectively.
From the values of model parameters given in Table 1, the calculated value of S is 0.14 and less than 1.63.The rupture velocity increases with rupture length and can be higher than the S-wave velocity or even approach the P-wave velocity.From Fig. 5b, the break tip is 1.8 km and time is 0.5 s, so the rupture velocity is 3.6 km s -1 , which is higher than the S-wave velocity.From Fig. 5c, the break tip is 3.7 km, so the rupture velocity is 3.7 km s -1 .Results suggest an increase in the rupture velocity with the rupture length.
The time-dependent displacements for the three faults are shown in Figs.6a -f.The maximum slip on the M-fault occurs at the center of the initial crack (Fig. 6a) and increases with time.From Figs. 6e and f, we can see that when the displacement of M-fault does not change significantly, the fault is almost stable and thus the slip ceases.The two sub-faults start to rupture about at t = 1.37 s.However, their displacements are not zero before they break because the waves generated from the M-fault passed through the two sub-faults.This can be seen from the panels on the left and right in Fig. 6c.The displacements on the two sub-faults are basically the same, but the patterns are asymmetric.Small differences are observed at t = 2.5 s.The differences are attributed to different amounts of later reflections from the free surface at different depths.The variations in the timedependent shear stresses on the three faults are shown in Figs.7a -f.Since the variations in normal stresses are very small, they are not displayed in Fig. 7.At t = 0, on the Mfault, the stress drops occur only at the site where the initial crack occurs and does not change at others.At t = 0.5 ~ 1.5 s, the maximum shear stress appears at the boundaries between the cracked and un-cracked regions, as displayed by blue in Figs.7b -d.These are crack tips which are ready to induce further cracks.At t = 2 and 2.5 s, the displacements within the M-fault are greater than the critical displacement (D c ) (see Fig. 3).In other words, when the slip distance is equal to D c , the shear stress is equal to the dynamic frictional strength, and thus, does not decrease further (Figs.7e and f).At t = 0 ~ 1 s, the shear stresses on the S-and N-faults change very little (see the left-and right-hand panels of Figs.7a -c).At t = 1.5 s (Fig. 7d), significant increases in the shear stresses caused by the ruptures on the M-fault are observed in the left-and right-hand panels of Fig. 7d and trigger the failure of the two sub-faults.At t = 2 to 2.5 s, the ruptures occur on the two sub-faults with a similar propagation pattern as happened on the M-fault.

DISCUSSION
To study fault interaction, it is very important to locate the minimum distance from which the main fault is capable of triggering sub-faults.This distance may depend on the mainshock magnitude, the geometry of the main fault and    where there are larger stress changes.The TRSR can be clearly identified based on the M-fault model.Based upon rupture-time contours and the initial triggered rupture area in Fig. 5f, we believe that if the sub-fault can completely rupture, most of the sub-fault must be within the area where the ruptures of the main fault stopped first.It is found that such an area is in between 0 and 1.73 km (measured from the western edge of the M-fault) on the N-fault and in between 1.73 and 3.46 km on the S-fault.But, the TRSR would be distinct for different models and model parameters.This means that when a sub-fault or any part of the two sub-faults in a to-be-triggered region, it could be triggered and would completely rupture.When the two sub-faults in Fig. 8a are shifted somewhat to the western portion of the M-fault, their eastern edges are still on the western edge of the Mfault.Hence, the two sub-faults are, respectively, inside and outside of the to-be-triggered regions.For this case, the Nfault ruptures completely, and the equal-time contours of ruptures become narrower at the place near -1 km (in the X-direction).However, it is not easy to interpret this phenomenon due to a complicated pattern of stress propagation.On the S-fault, only a local portion breaks and ruptures do not propagate over the entire fault plane.The incomplete rupture on the S-fault is due to low stresses exerted on it caused by the rupture on the M-fault.
In order to shift the two sub-faults from the M-fault as mentioned above, we define the distance of 0.865 km to be a length unit.First of all, the two sub-faults are shifted 3.46 km (or 4 length units), their respective eastern edges are on the western edge of the M-fault, move the two subfaults eastward by step in 1 unit length and calculating the triggering rupture until their respective western edges are on the eastern edge of the M-fault, as shown in Figs.8a -i.Figures . 8a -g show that the N-fault is located in the north part of TRSR, and can be triggered with a complete rupture.But, in Figs.8h and i, the N-fault is outside the TRSR, it will be triggered with an incomplete rupture.On the other hand, when the S-fault is located in the south part of TRSR, it can be triggered with a complete rupture, as shown in Figs.8c -i.When the S-fault is outside the TRSR, it will be triggered with an incomplete rupture, as shown in Figs.8a and b.
Obviously, simulations are carried to explore the triggering effects resulting from the related positions, along the east-west direction, which are perpendicular to the strike direction of the M-fault to the sub-faults.Meanwhile, the simulations were also conducted to study the effects due to the change of the distances between the M-fault and the subfaults along the north-south direction which is parallel to the strike.Results are shown in Fig. 9.In Figs.9a -e (i.e., the left and right panels for the S-and N-faults, respectively) we can see the modes of complete, incomplete, and non-triggered ruptures for five different distances, i.e., 1.2, 1.6, 2.0, 2.4, and 2.8 km, away from the M-fault in the north-south direction.The rupture pattern of the M-fault is the same as before and is not shown here again.Simulation results clearly show that the two sub-faults can be triggered with complete ruptures only when the separation distances are shorter than 2.0 km (Figs.9a -c) and with partial ruptures when the separation distances are in between 2.0 and 2.8 km (Fig. 9d).On the other hand, the sub-faults cannot be triggered when the separation distance is longer than 2.8 km (Fig. 9e).Of course, the separation distance in controlling triggering is not fixed at the value of 2.8 km, and depends on the models with different directions and magnitudes of regional stresses.A key issue that cannot be ignored is the direction of the applied regional stress.According to our simulations, the two sub-faults of the model can be triggered only when the angle φ lies between 18° -42° as mentioned by Wang (2005).We found that the TRSR is no longer located on either the northern or southern side of the M-fault when φ < 18° or φ > 42°.This means that the changes of shear stresses in the east-west component caused by the break of the M-fault are too small to trigger the two sub-faults.

CONCLUSIONS
In this study, we first develop and test a composite finite difference method to study the triggering and propagation of ruptures in a complex fault system with a main thrust and two strike-slip-type sub-faults which can be observed in central Taiwan.Two different grid sizes of the finite differences are performed for two different coordinates which are connected to each other by using a 2-D interpolation and coordinate transformation.This makes simultaneous computations for different faults geometrically possible thereby facilitating the study of fault interactions in such a complicated fault system.
The direction of applied regional stress does, of course, play a very important role in fault interactions, but this is, to a certain extent, limited.With a fixed direction of regional stress (φ = 26.57° in this study), the present simulation results provide solid unambiguous evidence that the strike-slip-type sub-faults can be triggered and then rupture completely by the failure of a nearby thrust fault under two conditions.The first condition is that the strike-slip-type sub-fault must be located inside the range of a source region (TRSR) for to-be-triggered by the ruptures of the main thrust.The second condition is that the strike-slip-type subfault is not too far away from the main thrust.When the two conditions do not exist on the strike-slip-type sub-faults, the rupture-time contours are narrow or the sub-faults cannot completely rupture even though they are triggered.In the present study, we simulate the rupture and triggering processes of a complex fault system composed by three faults.But, the real fault system is much more complicated than the present one.Hence, more studies will be conducted on the basis of numerical simulations to investigate the key factors, including the fault size, geometric patterns and ini-tial regional stresses, in influencing the triggering distance in the near future.
to the Japan Geological Survey.The authors thank two reviewers for giving significant suggestions and comments to substantially improve the article.This study was partly sponsored by Chia-Nan University of Pharmacy and Science.From a similar derivation, the displacement at the grid point on the footwall, i.e., ( , , , ) , is: Fig. 1.(a) The fault system is composed of a main thrust (the M-fault) and two strike-slip faults (the N-and S-faults).The faults are 4 km below the free surface.The main thrust with a strike in the Y-direction and a dip angle of 30° in the +X-direction is on the F-coordinate ( Xl Yl Zl ), while the other two sub strike-slip faults are on the O-coordinate (XYZ).(b) Plane view of the three faults is given in (a), with the two horizontal compress stresses, i.e., XX v and YY v , and the vertical stress ZZ v .

Fig. 3 .
Fig. 3.The slip-weakening fracture criteria for shear stress, v , 0 v is the initial shear stress; s v is the maximum static frictional stress; d v is the dynamic frictional stress; and D c is the critical displacement.When v > s v , the slip starts.When the slip is less than D c , v decreases linearly with increasing slip.When the slip reaches D c , v = d v .

v
of normal stress on the sub-fault at time t; shr 0 is the initial shear stress on the sub-fault; and ( ) t shr v is the change of shear stress on the sub-fault at time t.

Fig. 4 .
Fig. 4. The relation between φ and s D on the two strike-slip sub-fault planes at t = 0 s.

Fig. 6 .
Fig. 6. Figure shows the time-dependent displacements of the three faults.There is the same description as mentioned in Fig. 5.The slip starts at the center of the M-fault (a) and reaches a steady value after the rupture is completed.The asymmetric displacement patterns are caused by upward (N-fault) and downward (S-fault) rupture propagations.

Fig. 7 .
Fig.7. Figure shows the time-dependent displacements of the three faults.There is the same description as mentioned in Fig.5.At t = 0, stress drop is only observed at the site where the initial crack occurs, while all other regions maintain regional stress.The stress drop extends outward with time, as shown in a to f.The maximum shear stress appears at the boundaries of the cracked and un-cracked regions, as indicated by the blue shading in b to d; these represent crack tips and are ready to induce further cracks.A similar phenomenon is observed for the N-fault (right side) and S-fault (left side).
. Here, we will change the positions of the two sub-faults, as shown in Fig.8, to investigate the possible effects of the spatial structure of the fault system triggering the two sub-faults.The northern and southern ends in each panel of Fig.8marked by two slashed rectangles represent the triggered ruptured source regions (denoted by TRSR)

Fig. 8 .
Fig. 8.The relationship between triggering with changes of the locations of the two sub-faults in the east-west direction.The two slashed rectangles in each figure represent the triggering regions and this indicates that triggering for completely rupture is only possible when the sub-faults lie in this range (see text).

Fig. 9 .
Fig. 9.The relationship between triggering with changes of the two sub-faults in the north-south direction.The distances between the M-fault and the sub-faults are: (a) for 1.2 km; (b) for 1.6 km; (c) for 2.0 km; (d) for 2.4 km; and (e) for 2.8 km.The complete ruptures can be observed only for distances less than 2.0 km.

Table 1 .
The parameters and their values used in numerical simulations.