Estimation of the Wenchuan Earthquake Rupture Sequence Utilizing Teleseismic Records and Coseismic Displacements

For the 12 May 2008 Mw 7.9 Wenchuan earthquake, two imbricate faults, Beichuan fault and Pengguan fault, have ruptured simultaneously. Special attention should be paid to the point of 40 km northeast of the epicenter, in which the Xiaoyudong fault intersects the above two faults, creating a complex fault structure. Surface rupture data from field surveys and previous research of dynamics studies indicate that an important transformation may take place at the intersection. But, few studies about inversion of source rupture process have focused on this issue. We establish a multiple-segment, variable-slip, finite-fault model to reproduce the rupture process and distinguish rupture sequence. Based on the nonnegative least square method and multiple-time-window approach, the spatial and temporal distribution of slip for three rupture sequences are exhibited, using teleseismic records and coseismic displacements.*e conformity between synthetic and observed teleseismic records as well as the slip value of the shallowest subfaults and the coseismic displacements is utilized to calibrate the model. *e results are as follows: (1) *e teleseismic records inversion alone could not distinguish different rupture sequences. However, in order to make the slip of the Hongkou and Yingxiu area coincide with the field investigation, only the Beichuan fault has a bilateral rupture on the point of intersection of Xiaoyudong fault. So the possible rupture sequence is that the earthquake started at the low dip angle part of southern Beichuan fault, and then it propagated to the Pengguan fault, which caused the rupture of Xiaoyudong fault. *en the southern part of Beichuan fault with high dip angle is triggered by the Xiaoyudong fault. (2) *e coseismic displacements constraint can control the slip of subfaults near the surface and has little impact on the deeper subfaults. (3)*emaximum slip on the fault is located near the Yingxiu and Beichuan area; moreover, the slip is mainly distributed at the shallow region rather than at the deep, which led to serious disasters. Meanwhile, majority of the aftershocks occur in the periphery of large slip.


Introduction
e May 12th 2008 M w 7.9 Wenchuan earthquake which occurred in Wenchuan, China, caused a large number of casualties and very serious engineering damage.Field surveys showed that the earthquake ruptured the middle segment of the Longmenshan fault zone; meanwhile, two rupture spaces were approximately parallel with the northeastern Beichuan-Yingxiu fault (abbreviated as BCF) with 240 km surface rupture and the Guanxian-Jiangyou fault (abbreviated as PGF) with 70 km surface rupture [1][2][3][4].e BCF was divided into two sections with different movements at the point of Gaochuan [5][6][7][8][9].e southern part of Gaochuan (southwestern BCF), the slip on the fault, is mainly thrust, while the northern part of Gaochuan (northeastern BCF) slip type changes into strike.e research of the focal mechanism of aftershocks showed that aftershocks of the northern segment of the BCF (northern part of Gaochuan) were characterized by high dip angles [10].e width of aftershock distribution in the northern segment of BCF was obviously narrower than that in the southern section [11], which indicated the dip angle of the northern segment to be larger than that of the southern section.
It should be noted that the Xiaoyudong fault (abbreviated as XYDF) with strike about 325 °and length of 6 km [12], which was a small tear fault almost perpendicular to the PGF and BCF, located at about 40 km northeast of the epicenter.
Near the intersection with the XYDF, surface rupture of PGF and BCF showed obvious dislocation and discontinuity [1,13] (Figure 1).Surface rupture of PGF and BCF was 5 km and 7 km apart near the southern part of XYDF, while near the northern part of XYDF was 7 km and 15 km apart.Aftershocks along the XYDF demonstrated a leftlateral and strike-slip; nevertheless, a thrust slip of BCF at south of Gaochuan, which indicated that an important transformation may have taken place at the intersection.According to the result of Wang and Liu-Zeng [14], the XYDF played a positive role in linking coseismic slip transfer between the BCF and PGF, which used a quasistatic stress analysis.is raises the possibility that rupture of the XYDF caused bilateral rupture of the BCF.Inversion of the source rupture process also indicated that the BCF was initiated at the surface near the intersection of the XYDF [15,16] although this view di ered from most previous studies that have dealt with this rupture process [3][4][5][6][7][8][9][17][18][19][20][21][22][23].ese research studies presumed a unilateral rupture.In other words, the rupture sequence of BCF, PGF, and XYDF has been ignored, which is a valuable research issue.Furuya et al. [20] used crustal deformation data sets to develop a fault source model of Wenchuan earthquake, which considered the XYDF.But the rupture did not involve the time dimension, so the rupture sequence could not be recognized.
In this paper, we present a three-dimensional (3-D) fault model combined with an investigation of the aftershock distribution, seismic re ection pro le, and surface rupture data.Based on teleseismic records as well as coseismic displacement data, we invert spatiotemporal variation of the rupture process and analyze possible rupture sequences between BCF, PGF, and XYDF.Here, two standards of the judgments are utilized to calibrate the model.One is the degree of conformity of synthetic and observed teleseismic records, and the other is the slip value of the near-surface fault and the coseismic displacements.

Fault Rupture Model
Study on seismogenic structure of the Wenchuan earthquake showed that from south to north, the dip angles of the Beichuan fault changed signi cantly, the dip increased from deep to shallow in the southwestern BCF, and the northeastern section of BCF had a greater dip angle than the southwestern [1,6,[24][25][26][27].Based on the above features, a 3-D fault model for the Wenchuan earthquake was established (Figure 2).We parameterize this complex fault using six planes, one segment for PGF and ve for the southwestern and northeastern BCF.For the southwestern BCF, from deep to shallow, the faults are labeled as BCF4, BCF3, BCF2, and BCF1, with dip angle 20 °, 33 °, 50 °, and 65 °.Northern BCF is BCF5 with dip angle 60 ° [21], and dip angle of PGF is 30 °.
BCF4 intersects BCF3 at a depth of 16.6 km, BCF2 intersects the PGF at a depth of 10 km, and BCF2 intersects BCF1 at a depth of 5.4 km. e length of aftershock distribution is obviously larger than the surface rupture and in the north and south ends exceed 50 km [6]. is suggested that there may be slipping without surface rupture at both ends.In order to nd the possible rupture area in the deep region, the length of the fault is larger than the surface rupture.e PGF and southwest BCF are 132 km, northeast of BCF is 180 km, and the total length of BCF is 312 km. e strike of PGF and BCF is 224 °according to the surface rupture direction [1].e hypocenter is 30.986°N, 103.364 °E, according to United States Geological Survey.e initial rupture point is on the BCF3 and the depth is 14 km (China earthquake networks center).e subfaults are 5 km long and 3 km wide, and total subfaults are 854 (Table 1).

Data
Two di erent data sets, teleseismic body waves and surface o sets, are used.e teleseismic waveform data are from the Incorporated Research Institutions for Seismology (IRIS) and are deconvolved from the instrument response, integrated to obtain displacements.We use the distinct P-wave displacement observed at 36 stations with epicentral distances between 30 °and 90 °(Table 2).e locations of these 36 stations are shown in Figure 3, and the uniform azimuthal coverage can limit slip on the fault.e records are band pass ltered with a zero-phase third-order Butterworth lter in the range between 0.02 and 0.5 Hz and resembled with a time step of 0.2 s. e PP-waves are not used.
Another important data set we used to constrain the slip is the coseismic displacements.Field investigation shows  Advances in Civil Engineering that surface breaks extend for about 240 km from YX to QC on BCF and about 70 km from the intersection of XYDF and PGF to HW on PGF, as the red curve shown in Figure 1.
According to the results of Xu et al. [1,3,4], the surface rupture showed that the slip on the southwestern BCF was dominated by thrust slip; the northeastern segment is mainly the strike slip.e surface rupture of PGF is concentrated in the region between the intersection of PGF and XYDF and HW, with mainly dip slip.Since the length of the subfault of our model is 5 km, more than one observed surface offset may exist in one subfault, we averaged the measured surface slip along the 5 km length of each subfault and assigned that value to the shallowest subfault element.Figure 4 shows the total slip on BCF and PGF according to the results of Xu et al. [1,3,4].

Inversion Method
After the earthquake, the slip of the fault plane generates the seismic waves which are recorded by the stations.In order to obtain the slip values of the fault plane, the inversion method is usually utilized to achieve this goal.For the inversion, the input is the observed records.Meanwhile, the output is the slip vectors of the fault plane.We apply the nonnegative least-squares inversion method to calculate the slip on the fault according to Hartzell and Heaton [29].e fault plane is divided into a series of subfaults with the same size.e multitime window method can simulate heterogeneity of amplitude, duration, and direction for slip in the fault, as well as the rupture velocity.Every subfault can slip in any time windows, and then a complicated rise time function can be constructed.We assume a constant rupture velocity; if the slip occurs in subsequent time windows, it is interpreted as a lower rupture velocity, otherwise equal to the set value.
Here, the slip is divided into two directions, dip and strike, and the sum of these two values gives total slip and direction.In this approach, the inversion function is formed by the synthetic waveforms (matrix of A), solution vector of slip on each subfault (matrix of x), and the observed data vector (matrix of b): e function can be solved by linear least-squares method, but the result is unstable because matrix of A is an ill matrix.is leads to a phenomenon that the matrix of A or b has a small change, and the results of matrix x show greater variations.
e problem can be solved by adding damped and constrained to the above function.
Here C −1 d represents a priori data covariance matrix which is used to set each record in the inversion accounts for the same weight.e matrix of M for minimizing moments is obtained by letting x i � 0. e smoothing constraints matrix of S is used to minimize the slip between adjacent Advances in Civil Engineering subfaults as well as adjacent windows for one subfault along strike and dip direction.e slip on the peripheral area with small values is achieved by the Matrix of B.
e e ect of Matrix F is letting the sum of all time windows for surface subfaults to equal coseismic displacements (matrix of D) at corresponding position [28,30].Weights λ 1 , λ 2 , λ 3 , and λ 4 control the trade-o between satisfying these constraints and tting the data, and the best results are obtained using the minimum value of waveform mis t, as follows: In this expression, x i and y i denote observed and synthetic waveform data.e Green function is calculated using the ray theory method with the velocity structure model Ak135 and the lter band and sampling interval are the same as records.For Wenchuan earthquake, ve isosceles triangle time windows are used; the duration is 2 s with total rise time 10 s and adjacent subfault with a delay of the same length.
e dip and strike slip are 90 °and 180 °according to the rake angle of this earthquake.

Rupture Sequence
According to the fault model in this paper and previous corresponding research results [14-16, 31, 32], in order to analyze the rupture sequences among BCF, PGF, and XYDF, three possible cases are used (Figure 5).In Case 1, the rupture was initiated at the low dip angle part of BCF3, and then it propagated to deep area BCF4 and to shallow area BCF2 and PGF.Rupture front propagated through BCF2 and caused the rupture of BCF1 (process 1).In Case 2, the rupture was initiated at the low dip angle part of BCF3, and then it propagated to deep area BCF4 and to shallow area PGF (process 1).e XYDF was triggered by the PGF, and the BCF1 was triggered by the XYDF.
en BCF1 and BCF2 produced bilateral rupture (process 2).In Case 3, the rupture was initiated at the low dip angle part of BCF3, then it propagated to deep area BCF4 and to shallow area BCF2 and BCF1, successively (process 1).In this sequence, the XYDF was triggered by the BCF1, and the PGF was triggered by the XYDF.en PGF produced bilateral rupture (process 2).

Rupture Velocity Sensitivity.
e rupture velocity of Wenchuan earthquake was quite di erent with the range 2.5∼3.4 km/s [5,9,16,17,33,34].To search the optimal value, the maximum rupture velocity is 2.0∼3.4 km/s, about 0.6∼1 of the shear-wave velocity in the source region with interval 0.2 km/s, and the rupture velocity has 8 values.e teleseismic waveforms are used to study the rupture history for the three rupture sequences.e rupture velocity sensitivity for rupture 1, 2, and 3 is shown in Figure 6.
e result shows that the change of residual for the three rupture sequences is similar.e residual has little change with velocity less than 3.0 km/s and obviously increases with velocity greater than 3.0 km/s.e residuals of 2.2 km/s and 3.0 km/s are smaller.e slip of the area near Nanba on BCF5 contributes most to the end of synthetic waveforms with velocity 2.2 km/s; therefore, the slip near Nanba is very small, and the value is inconsistent with the observed data.By contrast, the slip near Nanba contributes to the synthetic records at about 60 s with large amplitude when the velocity is 3.0 km/s, so the slip is about 2.0 m and in good agreement with the observed value.Furthermore, the Wenchuan earthquake has a large rupture scale and complicated rupture process, and the range of rupture velocity may be larger.In order to simulate the rupture process, the velocity  Advances in Civil Engineering 3.0 km/s is set as the highest value.e residuals with velocity 3.0 km/s of rupture sequence 1, 2, and 3 are 0.271, 0.268, and 0.267, respectively, very close to each other.So, the teleseismic data alone are not capable of choosing the best rupture model without additional constraints.e same conclusion was obtained by Hartzell et al. [16].

Coseismic Displacements Constraints.
Figure 7 shows the comparison between the surface o sets and inverted slip of the near-surface fault with velocity 3.0 km/s for three rupture sequences.On PGF, the inverted slip in the south of BL is quite small for sequences 1 and 2, which agree well with the observed value.However, the inverted slip value of sequence 3 is about 2 m and the maximum is 3.8 m, which is seriously overestimating the true value.e reason is that the rupture process of sequence 3 is di erent from 1 and 2 on PGF.For sequence 3, the PGF has a bilateral rupture from the point of intersection with the XYDF.en, the south surface area of PGF has the same initial rupture time with the north area.e distance between the two regions is about 50 km, and the far eld Green function is similar.For the teleseismic inversion, when the Green function, initial rupture time, and the travel time of seismic phase are close for di erent area, the slip cannot be distinguished.So the slip in the northern part of PGF is transferred to the southern part.e slip inversions using teleseismic data only are a ected by a trade-o between rupture timing and slip location [16].e slip near the surface area of southern part of PGF contributes to the synthetic waveforms at about 30 s, so the slip is larger.On the contrary, the slip of sequences 1 and 2 for the same area contributes to the synthetic records at the beginning of wave packet, so the slip is very small.ere is a quite large surface o set between YX and HK areas on BCF. e result of sequence 2 is close to the observed value, but the inverted slip of sequences 1 and 3 is massively underestimated.At the point near HK with maximum surface o set, the value of sequence 2 is about 3/4 of the observed data.Nevertheless, the result of sequences 1 and 3 is only 1/7 of the observed value.e inversion results of the three sequences are similar, except the above two regions.rough analysis of the di erence between inversion results and coseismic displacements, it is shown that the surface area of the southern PGF of sequence 3 has a large slip value and the area from HK to YX on BCF of sequences 1 and 3 has a small slip value, which are contrary to the observation.So sequences 1 and 3 are not reasonable.For these two areas, the result of sequence 2 is in good agreement with the observation.From the surface o sets, it is considered that the rupture sequence 2 is in accord with the true rupture process of Wenchuan earthquake.For the rupture process of earthquake, the situation which the propagation is to the reverse direction of rupture is infrequent.However, the same phenomenon was proved to exist.e 1984 M6.2 Morgan Hill California earthquake, 1992 M w 7.3 Landers earthquake, and 2010 M w 7.2 EI Mayor-Cucapuh earthquake had the same rupture mode [35][36][37].is phenomenon is coincident with the numerical simulations of the rupture process [38].
e complex structure of seismogenic fault leads to such a rich rupture process.

Spatiotemporal Variation of the Rupture Process of the Wenchuan Earthquake.
Comparison between observed teleseismic records and synthetic records of three sequences, as well as the slip distribution of rupture velocity 3.0 km/s, is shown in Figure 8.For the three rupture sequences, the slip is concentrated at 5 areas.e rst, near the initial rupture point, the area is small but the thrust slip is large.e second, under LMS on BCF3 and BCF4, the slip block is large (A 1 ).e slip is mainly the thrust slip on BCF3 and the right lateral strike slip on BCF4.e third, from YJS to GC on BCF1 and BCF2, there is a large slip area (A 2 , A 2-1 ).At the surface area, it is mainly the right strike slip, while on the bottom area is dominated by the thrust slip.e fourth, on BCF5, the slip at surface area near BC is the large thrust slip (A 4 ).e fth, on BCF5 between NB and QC, the slip is thrust and strike (A 5 ).On the other side, slip distribution is obviously di erent for three sequences at some areas.e rst, on the high dip angle part between HK and YX of the southwestern BCF, sequences 1 and 3 have no slip.Nevertheless, for rupture sequence 2, a large slip patch is located in the region dominated by mainly the thrust slip.e second, on PGF, slip distribution is similar for sequences 1 and 2, a large slip patch is located at the bottom of the fault near the middle area, but nearly no slips exist on the southwestern PGF.In contrast, for sequence 3, the large slip patch shifts to the southwestern PGF.
e synthetic records of three rupture sequences show good agreement with the observed ones.Because of the close spatial distribution of the faults and the similar steep takeo angles, the synthetic teleseismic records are almost same for the three rupture sequences.is leads to the fact that teleseismic data alone cannot distinguish between di erent rupture sequences, but by combining the result of fareld record inversion and the surface o sets, the optimal rupture scenarios can be distinguished.

6
Advances in Civil Engineering In order to limit the slip of near-surface subfaults, the surface o sets are involved in inversion.
e result is shown in Figure 9, the Figure 9(a) shows the slip distribution, Figure 9(b) shows the moment rate function for the complete rupture, Figure 9(c) shows the rupture process snapshots with 4 s interval, and Figure 9(d) shows the comparison of sevenmonth aftershock distribution and the projection of the spatialslip distribution.Table 4 comprises the parameters of asperities.
e total seismic moment is 0.950 × 10 21 N•m, which is close to the result of USGS and Harvard CMT solution.e result shows that the slip along the southwestern BCF is near the initial rupture point and encompassing the area between YX and LMS (asperity A 2-2 ), YJS and GC (asperity A 2-1 ), and the area under the LMS (asperity A 1 ).e corresponding moment is 0.387 × 10 21 N•m, accounting for 41% of total seismic moment.e slip on PGF is located near Bailu (asperity A 3 ), and the corresponding moment is 0.115 × 10 21 N•m, 12% of the total value.Generally, the slip on the southwestern BCF and PGF is dominated by the thrust.For BCF5, the slip is located near the surface of BC (asperity A 4 ), as well as in the area between NB and QC (asperity A 5 ).e moment of BCF5 is 0.446 × 10 21 N•m, 47% of total moment.
Our inversion result shows that the slip is mainly distributed at the depth above 15 km.At the fault plane, there are 6 asperities, which indicates that the earthquake is composed of at least 6 subevents.
e slip is mainly distributed on the BCF, indicating that the BCF is the main rupture fault.According to the moment rate function and rupture process snapshot, the whole rupture process about 100 s is divided into 4 stages.0∼12 s is the rst stage, with 3% of the total seismic moment, corresponding to the slipping near the initial rupture point.e second stage, from 12 to 40 s, is the process of concentrating release of energy, corresponding to the rupture of southern segment of BCF and PGF, releasing 50% of the total seismic moment.As shown in Figure 9(c), at about 12 s, the rupture front arrives at the area under LMS, and then asperities A 1 and A 3 begin to rupture.At 16 s, asperities A 2-1 and A 2-2 begin to produce a large slip; meanwhile, there is a bilateral rupture on the BCF from the point of intersection with the XYDF.e third stage, from 40 to 60 s, corresponding to the rupture near BC, released 18% of the total moment.60∼88 s is the last stage, which releases 25% of the total moment, corresponding to the rupture of the region south of BC.
e outstanding peculiarity of the moment rate function is that one-half of the energy is released in the second stage which lasted 28 s, a short duration for an earthquake of that size.In such a short period of time, a great deal of energy is released, causing serious damage.
rough the analysis of aftershock distribution and slip region, some useful results are obtained.One is that majority of the aftershocks occurred in the periphery of large slip.It shows that most of the stress was released in the area of asperities so that no aftershocks occurred after the main shock.e pattern of aftershock location in relation to the areas of large slip is coincident with the 1999 M w 7.6 Chi-Chi earthquake and 2015 Nepal M w 7.9 earthquake [39,40].
ere are seven aftershocks with magnitude greater than 6.0 around the fault.It is remarkable that all the aftershocks larger than 6.0 occurred near the large slip region.Two aftershocks were located near the epicenter, and three were located at the periphery of the northernmost region.e other is that the aftershocks tend to be distributed at the area of the projection of the deep part of the fault.From the research of Jia et al. [27], the most dense aftershocks were concentrated from 14 to 16 km.e result of Tong et al. [8] and Qi et al. [22] showed the same situation.It is suggested that for the Wenchuan earthquake, the slip is mainly distributed at the shallow region, while the aftershocks mainly occurred in deep regions, which led to fateful disasters.

Conclusion
e rupture process of Wenchuan earthquake is complicated including an obvious change in slip direction, which suggests that several faults participated in the rupture.Especially at the area of 40 km northeast of the epicenter, the XYDF intersects the BCF and PGF, forming complex fault geometry.e surface ruptures and results of dynamic studies indicated that an important transformation might take place at the intersection.Unfortunately, few research studies about inversion of source rupture process have focused on this      e red star is the epicenter, and the brown stars are the epicenter of seven aftershocks with magnitude greater than 6.0 [10].

Figure 4 :
Figure 4: e slip values of shallowest subfaults of BCF and PGF according to the results of Xu et al. [1, 3, 4].

Figure 5 :
Figure 5: ree rupture sequences for the Wenchuan earthquake: (a), (b), and (c) refer to rupture sequences 1, 2, and 3, while the circles denote the earliest rupture subfaults, and the arrows denote rupture directions.

Figure 6 :
Figure 6: Rupture velocity sensitivity for rupture sequence 1, 2, and 3. e velocity is tested from 2.0 to 3.4 km/s.e preferred model has a value of 3.0 km/s.

Figure 7 :
Figure 7: Comparison between inverted slip of the near-surface fault and coseismic displacements from eld survey.Bar graph is the surface o sets.

Figure 8 :
Figure 8: Comparison between observed teleseismic records (black) and synthetic records of sequence 1 (red), 2 (blue), and 3 (brown).e slip distribution of rupture velocity 3.0 km/s of rupture sequences 1, 2, and 3. Point A and B are the surface intersection between the XYDF and PGF as well as the southern BCF.

Figure 9 :
Figure 9: (a) e slip distribution with coseismic displacements constraints of rupture process 2. (b) Moment rate function.(c) Rupture process snapshots of Wenchuan earthquake with 4 s interval.(d) e comparison of seven-month aftershock distribution (solid circles) and the projection of the spatial-slip distribution of BCF. e information of aftershocks is from the China earthquake network center (CENC).ered star is the epicenter, and the brown stars are the epicenter of seven aftershocks with magnitude greater than 6.0[10].

Table 3
[28]orate the surface offset data to constrain slip on the shallowest area because we know where and how much the surface slipped[28].

Table 2 :
Locations and epicentral distances of 36 teleseismic stations.

Table 3 :
Slip values constraint of shallowest subfaults of BCF and PGF.

Table 4 :
e parameters of asperities.