Reconstruction of the Kinematics of Landslide and Debris Flow Through Numerical Modeling Supported by Multidisciplinary Data: The 2009 Siaolin, Taiwan Landslide

Chien-chih Chen1, Jia-Jyun Dong2, Chih-Yu Kuo3, Ruey-Der Hwang4, Ming-Hsu Li5 and Chyi-Tyi Lee2 1Grad. Inst. Geophys. & Dept. Earth Sciences, Nat’l Central Univ., Jhongli, Taoyuan 2Grad. Inst. Applied Geology, Nat’l Central Univ., Jhongli, Taoyuan 3Research Center for Applied Sciences, Academia Sinica, Nankang, Taipei 4Dept. Geology, Chinese Culture Univ., Taipei 5Grad. Inst. Hydrological & Oceanic Sciences, Nat’l Central Univ., Jhongli, Taoyuan Taiwan


Introduction
Typhoon Morakot struck southern Taiwan in the summer of 2009, causing the region's most severe flooding since the 1950s. In the early morning of August 9 (local time), a rainfalltriggered landslide and debris-flow extinguished the township of Siaolin Village, Kaohsiung. The interviews of many survived villagers could shed light on some parts of the story about the landslide catastrophe ). Scientifically, a simulation of the landslide/debris-flow can be used to examine its complex kinematic characteristics. A landslide/debris-flow simulation is conventionally performed using constrained digital terrain models, field geological mappings of channel cross-sections or laboratory measurements of slope material properties. Here, we employ seismological and near-surface magnetic data in a novel way to validate and reinforce our simulation of the catastrophic Siaolin debris avalanche. Siaolin is identified from a high-resolution aerial photograph taken before Typhoon Morakot (Figure 1). Changes in elevation before and after Typhoon Morakot ( Figure 1) were derived from two digital terrain models (DTMs) with a precision of five meters by the Agriculture and Forestry Aerial Survey Institute (AFASI) of Taiwan. Based on the DTMs before and after Typhoon Morakot, the major landslide body (cool colors in Figure 1) had an extent of 57 hectares, and the volume of the landslide was about 23.87 million m 3 . According to extensive field investigation, the sliding block was mainly composed of Pliocene shale and Quaternary colluvium. The sliding surface may have been located along the interface between the fresh and weathered Pliocene shale. The volumes of the deposits (warm colors in Figure 1) on the unnamed creek and the west bank of the Cishan River were estimated to have loose-measure volumes of 10.91 and 4.53 million m 3 , respectively. Most of the deposits laid down in the Cishan River, which had dammed the river for a short time after the debris avalanche, were washed away. The deposits that were washed away had an estimated volume of 10.83 million m 3 , assuming that the loose-measure volume of the sliding body yielded an expansion of 10%. Consequently, the total bulk of the natural landslide dam (i.e., the sum of the washed-away materials and the deposits left on the west bank of the Cishan River) was 15.36 million m 3 .

Frictional experiments of the sliding materials
Frictional properties of sliding surfaces at low to high slip rates are crucial for conducting numerical simulations of landsides/debris-flows. Recently, researchers have utilized the rotary-shear high-velocity friction apparatus to measure the frictional coefficients of the sliding surfaces of landside materials under common conditions present in landslide motion (Mizoguchi et al., 2007;Ferri et al., 2009;Miyamoto et al., 2009;Yano et al., 2009;Togo et al., 2009). The high-velocity friction apparatus and the experimental procedures we utilized in this study are described in papers by Mizoguchi et al. (2007) and Togo et al. (2009). Shale and colluvium in the Siaolin area were carefully sampled. The sampling locations are indicated in Figure 1. The tested samples were sustained under a normal stress of around 1 www.intechopen.com MPa, which accurately represents the in-situ normal stress for these sliding surfaces. For the colluvium samples, the initial water content was kept equal to the natural water content of the colluvium collected in the field (~20%). For the wet shale samples, we mixed 0.8 g of shale powder with 0.2 g of pure water, corresponding again to an initial water content of 20%. A wide velocity spectrum from 0.01 through 1.3 m/s, covering the range of likely sliding velocities, was adopted. Figure 2A shows the frictional behavior of the wet colluvium at different slip rates. Except at 0.01 m/s showing almost constant friction coefficient of 0.3, friction usually decreased with shear displacement after its peak was reached, and the residual friction fell in the range of 0.05 to 0.2. The 'V-shaped' friction curve, i.e., slip weakening then slip hardening, observed in some experiments at 0.7 and 1.3 m/s may be due to the fluid pressure buildup followed by the fluid escape from the slip zone. Figure 2B displays the frictional behavior for the wet shale samples at three different slip rates. All showed that the friction decreased with shear displacement from the peak value of 0.2 -0.3 to its residual value of 0.15. Again, the friction at 1.3 m/s showed the 'V-shaped' behavior at displacement of 30 -40 m. This indicates the dried shale gouge, after the release of fluid pressure, may have quite low frictional strength at such a high slip rate of 1.3 m/s. In both cases of colluvium and shale, the residual friction coefficients are mostly around 0.1 -0.2. The experimental results thus suggest that a frictional coefficient of around 0.2 is a reasonable estimate for the relatively high-velocity sliding regime before the moving block collapsed and transformed itself into a granular flow. Although the frictional coefficient of a granular flow is difficult to determine experimentally, an equivalent frictional coefficient of ~0.2 also falls within the range between 0.1 and 0.2 inferred from theoretical calculations (Okura et al., 2000).

Simulation of Siaolin landslide
The apparent friction coefficients determined in the laboratory tests are helpful in the applications of continuum landslide models because they can be directly associated with the friction force. In the present study, we use a two-dimensional shallow-water model extended from the paper of Tai and Kuo (2008) to simulate the Siaolin landslide. The effects of the general topography on the landslide flow are brought into account. The flow is simplified as an incompressible and inviscid fluid with a basal Coulomb friction force. The sole rheological parameter, the friction coefficient, is assumed to be a constant. Though with such simplifications, this type of models has been successively validated to be able to reproduce experimental and natural granular flows (Gray et al. 2003;Kuo et al. 2009). In addition to be constrained by the frictional experimental results, the Coulomb friction coefficient in the simulation is further referenced to a universal scaling law (Staron and Lajeunesse, 2009) and reconfirmed by an iterative optimization scheme. Inspecting extensive practical landslide data and discrete element simulations, the friction coefficient is found to follow a simple geometrical relation, the universal scaling law, which is inversely proportional to the one-third power of the landslide volume (Staron and Lajeunesse, 2009). Based on the friction angle 6 o found in the Tsaoling landslide with a volume 0.126 km 3 (Kuo et al., 2009), we estimate the friction angle for the Siaolin landslide (about 1/6 in volume) is about 11 o according to the geometrical scaling relation. This value is well constrained in the range suggested by the aforementioned frictional tests (assuming the friction coefficients asymptote to constant values). The initial volume is taken from the DTMs without a presumed volume dilation. The major slid mass on the east side of the Cishan River is www.intechopen.com released from the rest at t = 0 ( Figure 1). The computation domain 3,710×2,220 m 2 is discretized into a 10×10 m 2 grid mesh and the simulation time is 180 seconds. With the www.intechopen.com best-fit algorithm to minimize the difference between the simulation and the field measurement, it is confirmed that the optimal value of the friction coefficient is 11.5 o , with just a minute offset from the straightforward estimation. Fuller theoretical and numerical details of the simulation have been reported in a designated paper (Kuo et al., 2011). Two snapshots at t = 32 and t = 60 seconds are shown in Figures 3A and 3B, respectively. These times correspond to when the landslide flow reached the 590 Height (cf. Figure 1) and the west bank of the Cishan River. Upon hitting the 590 Height, the flow split into two sliding courses ( Figure 3A). The main stream was diverted to flow along the valley of the unnamed creek, and the second flow moved through the southern part of the 590 Height into the Siaolin Village. The simulation shows that when the flow hit the west bank of the Cishan River, a portion of the volume of debris overflowed from the guarding northeastern ridge into the Siaolin Village about 60 seconds after the initiation of the landslide ( Figure 3B). The velocity of the landslide's front reached a magnitude of about 50 m/s, and the duration of the landslide was about 110 seconds. The simulated deposit volume ( Figure 3C) in the Cishan River was about 21.1 million m 3 ; this formed a dam that was subsequently flushed away when it was breached about 30 minutes later. An animation of the whole process of the Siaolin landslide can be found in the online supplementary material.

Geophysical data validating the landslide simulation
The scenario replayed in the simulation agrees with two aspects of the data obtained from seismographs and a near-surface magnetic survey. First, several broadband TCWBSN (Taiwan Central Weather Bureau Seismic Network) seismographs recorded ground motions, especially Rayleigh waves, generated by the Siaolin landslide. We used the arrival times of the 25-s Rayleigh waves recorded at 8 stations to locate the original position of the landslide. A half-space model with 3.75 km/sec of S wave and an inversion algorithm (Stein and Wysession, 2003) were adopted for the determination of the landslide location. Estimated from the Rayleigh-wave calculation is the landslide location (Figure 4), together with its occurrence time of 6:16 in the moring. We note that at least three impact sub-events (marked by the arrows in Figures 4A and 4D) during the landslide can be clearly identified by the Rayleigh waves recorded, for example, at the distant Station YHNB. The seismogram recorded at Station SCZB (blue line in Figure  4A) also shows the visible sub-events during the landslide. Having the waveforms as recorded by TCWBSN, we simulated the waveforms of Rayleigh waves at a frequency band of 0.033-0.05 Hz. The focal mechanism of a horizontal thrust faulting, i.e. an equivalent single-force representation of shallow landslide, was used in the waveform modeling (Dahlen, 1993). The velocity structure is again a simple half-space model with 6.5 km/sec of P wave, 3.75 km/sec of S wave and 2.8 g/cm 3 of density. The program of synthetic waveform of surface wave was developed by Wang (1981). Also, to avoid the interferences from the velocity structure and the rupture directivity of source, we simulated the waveforms at the closer Station SCZB where is located to the south of the landslide. Three parameters in the forward modeling of the occurrence time, duration and moment for each sub-event were fitted to ensure that the synthetic waveform (red line in Figure 4A) matches the observed one. Three sub-events thus obtained in the modeling can fit well the waveform observed at Station SCZB indeed. Seismic moments for these three sub-events are 6.0, 1.4 and 3.0 10 16 Nm, respectively ( Figure 4B). This calculation shows that the second and third sub-events occurred 34 and 55 seconds, respectively, after the occurrence of the first sub-event. If we assume that the first sub-event represents the landslide's initiation, these two moments are quite coincident to the above-mentioned moments when the landslide flow reached the 590 Height and the channel of the Cishan River. After calculating the debris fluxes at two benchmarks near the 590 Height (AA' in Figure 3A) and the Cishan River (BB' in Figure 3B), we find that the large impacts to the 590 Height (blue line in Figure 4C) and the Cishan River (green line in Figure 4C) occurred around 40 and 60 seconds after the landslide's initiation, respectively. The waveforms recorded at Station SCZB are thus consistent with the kinematics resulting from the simulation of the Siaolin landslide event. In addition, the results of the near-surface magnetic survey ( Figure 5) provide another confirmation of the landslide simulation regarding the flow direction of the debris. Highresolution near-surface magnetic data covering the main township of Siaolin Village were obtained in October. The track lines of magnetic measurement were along the direction of NW to SE (with an interval of five meters). Portable proton precession magnetometers (Geometrics model G-856, with a sensitivity of 0.1 nT) were used in the magnetic surveys. The total-field magnetic data were collected at 1.22 meters and 1.82 meters above the ground. To ensure measurement reliability, the magnetic data were measured at least three times at each location within a few minutes. We refer readers to the paper of Doo et al. (2011) for the detailed information about magnetic survey and data analyses. The magnetic anomaly map (colors in Figure 5A) shows several anomalous patterns of magnetic dipole fields, and it can be used to derive the amplitude distribution (contours in Figure 5A) of the zeroth-order analytic signals of the magnetic sources (Hsu et al., 1998;Hsu, 2002;Doo et al., 2009). When superimposing the distributions of magnetic sources onto the aerial photograph taken before Typhoon Morakot, we can clearly see that all of the magnetic sources correspond well to the buildings located in the southern part of the township ( Figure 5B). Next, a sharp boundary of the destruction front of the landslide (thick, dashed lines in Figure 5) can be clearly identified by a comparison between the magnetic anomalies and the aerial photograph. Buildings located to the north of the boundary have been swept off their original sites by a massive impact. The landslide simulation of the Siaolin event suggests the existence of such a destruction front indeed, and its location coincides well with that boundary identified by the magnetic survey. Shown in Figure 5C is the simulated flowfield (arrows) of the landslide at t = 54 seconds in the corresponding region where the magnetic survey was conducted. We find that the impact of the landslide to the south, beyond the identified destruction front, is quite small. The final deposit heights (colors and black contours in Figure 5C) over the southern part of the village are less than one meter, suggesting that many of the buildings to the south may have been spared. The results of the filed magnetic survey thus strongly support the simulated flow pattern of the twodimensional shallow-water model. The magnetic results also draw the second phase of the disaster that buried the whole Siaolin -the fatal mud flow that followed the breach of the landslide dam. Note that within the area where the magnetic survey was conducted, the photo has been filtered to emphasize the buildings. (C) Simulated velocity (arrows) of debris 54 seconds after the landslide initiation and deposit height (colors and black contours) 74 seconds over the same area of (A), again superposed by the magnetic analytic signal (red contours). The thick red dashed line in each panel depicts the sharp destruction front of the debris suggested from the results of the magnetic survey.

Scenario of the breached landslide dam
While the landslide mostly destroyed the northern part of the Siaolin Village, the spared southern village obviously suffered different catastrophic scenario induced by the breach of the landslide dam. For a realistic simulation of the mud flow induced by the dam breach, the US National Weather Service BREACH model (Fread, 1991) was applied. This model generates a landslide dam break hydrograph, which in turn is used as input to drive the Federal Emergency Management Agency-approved FLO-2D model (O'Brien et al., 1993;O'Brien, 2006). The equivalent landslide dam in the model had an along-river width of 300 m and a crest elevation of 420 m. The equivalent cross-valley length of the crest was 213 m after compensating for the differences between the real valley topography and the simplified trapezoidal dam shape required by BREACH. The corresponding debris volume over the main channel of the Cishan River that failed is thus 5.2 million m 3 . The model geometry of the landslide dam has been clearly elucidated in the paper of Li et al. (2011), together with the detailed description about the simulation of the dam breach and mud flow. The initiation time of the dam failure was 6:43 a.m., and it only took about 8 minutes to reach a complete breach. The peak discharge rate of this massive breach was 94,280 m 3 /s. The average sediment concentration by volume was 0.38, which was estimated as the ratio of landslide dam volume to the total mud volume of the break hydrograph. Figure 6 depicts the deposited sediment depths over Siaolin as simulated by the FLO-2D model. Muddy floods following the breach covered the entirety of the village around 6:51 a.m., resulting in tragedy for the inhabitants of Siaolin. The deposited mud brought by the breaching dam is around 10 meters thick above the southern village, where was little ruined by the landslide. In terms of spatial distribution and depth, the overall pattern of the deposited sediment obtained from the elevation differences (contours in Figure 6) of the DTMs before and after Typhoon Morakot was well simulated by the model (color patches in Figure 6). Also note that the main river course was slightly diverted to the east, indicating that a significant scour occurred after the dam break deposition. Such a process of mixed scour and deposition is beyond the limitations of current FLO-2D modeling.

Concluding remark
Typhoon Morakot, an example of extreme precipitation, dropped more than 2,500 mm of rain over southern Taiwan within the three days following August 7, 2009. The National Disasters Prevention and Protection Commission in Taiwan reported 724 deaths caused by the consequent flooding. Of these victims, 474 were buried alive by the landslide and debrisflow event at Siaolin Village. This disaster had a severe socio-cultural impact on Taiwanese aboriginals of the Pingpu group. Reconstruction of the burying course due to the debris avalanche event will not only bring emotional satisfaction to the surviving villagers but will provide a fundamental understanding of the complex mass movement process (Densmore et al., 1997;Densmore and Hovius, 2000;Stark and Stark, 2001;Bachmann et al., 2004;Bruckl and Parotidis, 2005;Bonnet and Crave, 2006). Contrary to the fairly precise rupturing simulation of earthquakes, the complex mass movement of landslides/debris-flows is rarely reconstructed by numerical modeling (e.g. Kuo et al., 2009;Crosta et al., 2004;Lin et al., 2005). A quantitative description of the mass movement process of a landslide is thus a challenging task compared to the well-developed techniques associated with earthquake rupturing simulations. The accurate reconstruction of landslide kinematics is generally extremely difficult because of its complexity and of the lack of detailed data on mass movement processes. A landslide/debris-flow simulation is conventionally performed using digital terrain models solely, but one can raise the question whether the results of such simulations are realistic. We here employ interdisciplinary seismological and near-surface magnetic data to validate the results of a simulation of the catastrophic Siaolin, Taiwan landslide triggered by Typhoon Morakot in 2009. We have demonstrated that a realistic reconstruction of the catastrophic Siaolin event can be achieved by a simulation of the frictional continuum debris model validated with the support of geophysical (seismological and magnetic) data. In the first stage of the Siaolin case, a debris avalanche swept away half the village while in the second stage a debris-blocking dam was breached about 30 minutes after it was formed. Then muddy floods shortly covered the entirety of the village. This study shows that the pattern of the Siaolin event can be interpreted as a hybrid of mass movement processes of an avalanche, a landslide dam and a dam breach. Fig. 6. Simulated deposited sediment depths (color patches) around Siaolin Village due to the muddy flood of the breached landslide dam ( Figure 3C). Contour lines show the elevation differences between two versions of DTMs before and after Typhoon Morakot.

Acknowledgments
The scientific team of Project Morakot is grateful to the National Science Council (ROC) for its research support from grant NSC 98-2745-M-008-013. Also involved in this project are