The Impacts of River Channel Blockages Caused by Sliding Embankment Collapses during Earthquakes

In recent years, river blockages induced by debris flows have occurred more frequently because of the impacts of climate change and earthquakes, especially high-intensity earthquakes [1], which are often accompanied by multiple collapses and landslides [2,3] that change the geological and geomorphic conditions and increase the scale of loose materials in gullies, which are prone to large-scale debris flows when rainfall is heavy [4,5,6,7,8,9]. River blockages occur when the abundant loose materials of debris flows move from tributaries into main rivers within a short time, and the main rivers lack the capacity to carry such materials downstream. This often leads to upstream inundation and flooding outbursts induced by debris-flow dam failures. This process is often termed the “hazard chain” associated with debris flows [10,11,12,13,14,15,16]. In 2021, an embankment landslide in the upper reach of a river in Atami City, Shizuoka Prefecture, Japan, was associated with debris flows that caused extensive downstream damage [17]. Thus, the laws regulating residential land development and embankment seismic strengthening were revised, and a new embankment regulation law (the “Residential Land Development and Embankment Regulation Law”) was enacted [18]. It requires that disasters associated with embankment failures are to be minimized by the designation of regulated areas and assessment of the disaster risks [19]. Debris flows, historically termed “mountain tsunamis [20],” have occurred since the Edo period (1603–1868) and triggered landslides and river channel blockages during (for example) the 1847 Zenkoji earthquake in Nagano Prefecture. Landslide dams formed but were later breached, associated with downstream flooding and debris flow damage (for example, Figure 1) [21]. The new law focuses on the disaster risks associated with embankment landslides, such as river channel blockages, water accumulation, and flooding [19].

To date, many large-scale landslides and embankment collapses attributable to the effects of earthquakes have been reported [23,24,25]. A major disaster developed during the Hanshin-Awaji Earthquake [26]; residential embankments in the Nigawa Yurino area of Hyogo Prefecture, Japan, slid and blocked adjacent river channels (Figure 2 and Figure 3) [27,28,29]. There was no significant down-stream flooding because the flow of the blocked river in the Nigawa Yurino area was low and remedial measures such as sediment removal were rapidly implemented.

The historical risks associated with sliding embankment collapses have been evaluated [30,31]. The studies of Asada [32,33] and those of authors who analyzed the Kushiro-oki and Noto Peninsula earthquakes [34,35] highlighted the risks of embankment collapses followed by river channel blockages during earthquakes. The earthquake in the Noto Peninsula on 1 January, 2023, triggered the formation of many landslide dams [36]. The revision of the “Act on Promotion of Sediment Disaster Countermeasures in Sediment Disaster Hazard Areas” of May, 2011, mandated emergency countermeasures to address possible large-scale landslide dam formations; this enhanced social awareness [37]. The risk of debris flow and downstream flooding increases when the regular river flow is high or countermeasures are not quickly implemented.

All of the seismic forces in play, the ground structure, and the material properties (especially those associated with volcanic geology) affect embankment collapse during earthquakes. Most researchers have primarily focused on the direct disaster risks posed by such a collapse [38,39,40,41]; to our knowledge, the secondary risks have received little attention. In particular, the risk of flooding after a landslide dam is breached is significantly affected by the dam shape. Thus, the deposition morphologies of river channel blockages (landslide dams) associated with embankment collapses require more attention [42].

Several predictive equations have sought to define the shapes of landslide dam depositions based on small-scale water channel experiments [43]. When several thousand cubic meters of semi-dry sediment collapse into a stream, the landslide dam may be temporary. However, the reproducibilities of the predictive equations [43] have been questioned; analyses of historical landslide dam behaviors revealed a need for improvement. For example, the correlation coefficients of equations predicting land-slide dam deposition lengths ranged from only 0.16 to 44,45].

Several methods have been used to predict the movements of collapsed soil masses on slopes [46,47], and numerical analyses predicting the diffusion ranges of landslide soil masses (which were assumed to be incompressible fluids) have been published [48,49]. However, the predictions differ significantly because the values of various parameters (for example, those describing fluidity during movement) were arbitrarily chosen [49].

We thus developed simple equations that predicted the shape (the deposition height and length) of landslide dams that form after earthquakes and heavy rain in Japan; we referred to historical cases [44,45]. We earlier verified the reproducibilities of these equations during heavy rain [45] by reference to the large landslide dams that formed in the Kii Peninsula during Typhoon 12 (Talas) in 2011 [50]. However, the utilities of these equations when analyzing landslide dams that form after earthquake-induced embankment collapses have not been verified. Previous studies [45] developed more detailed predictive equations by reference to the geologies of collapsing mountainsides; they sought to improve the accuracies of deposition shape predictions after earthquakes. However, given the differences in the geologies of the collapsing embankments, application of the detailed predictive equations remains challenging.

Generally applicable, empirical, predictive equations derived via statistical analyses of international data have been used to predict outflow volumes when landslide dams that develop after river blockages breach [51]. However, these equations employ logarithmic relationships (log-log graphs) and are associated with predictive errors. Detailed individual examinations (numerical calculations) that refer to the specific conditions in play are essential when analyzing how the disaster risks are affected by the outflow volumes of breached landslide dams. Simple models predicting the effects of breaching have been proposed, but further verification is necessary [52,53].

In the context of Japanese landslide dams, models that combine one-dimensional outflow predictions with simple quasi-two-dimensional models have been developed to predict downstream flood damage when a dam overtops, and to facilitate emergency warnings [37,51]. However, these seek principally to predict outflows in river channels well downstream of the blockage point; the detailed effects in areas immediately downstream or upstream of the blockage point are not well-handled. Also, the models have not been verified in the real world; further work is necessary.

Here, we present an analytical model, the reproducibility of which, in terms of predicting dam flood outflows, has been verified by reference to landslide dam breach experiments performed in real-world streams [54,55]. We verified model reproducibility for steep, but not gentle, landslide dams, that is, those forming after blockages of river channels with steep or gentle slopes (Figure 2 and Figure 3). For the latter, the deposition height is low compared to the deposition length [55].

We here discuss the impact (disaster risk) of river channel blockages (landslide dams) triggered by the sliding collapse of embankments near rivers during earthquakes. We newly verify the reproducibility of simple equations that predict river channel blockage shapes [45,50] using examples from the Hanshin-Awaji Earthquake. Further, we confirm the reproducibility of our proposed two-dimensional model [55] for landslide dams with gentle slopes and we predict the dam outflow volumes and the disaster risks in the event of dam formation and breach in the Nigawa Yurino area.

To simply and accurately predict the deposition shapes of landslide dams created via large soil mass movements that block river channels, we now provide an overview of our earlier predictive formula, as shown in Figure 4 [45].

Generally, methods that predict the scale and deposition shape of a landslide dam include one or more physical models of mountain collapse and soil mass movement; both computational models [48,49] and statistical analyses have been used to analyze the relationships between historical deposition shapes and other factors. In our previous study [45], we used a statistical method to extract factors that correlated highly with the deposition height and length. The data (e.g., topographic conditions) were readily available; all were recorded during historical emergencies. We developed a predictive formula via the statistical processing of such factors. We first used the Quantification Theory Type III [56] of multivariate analysis to define trends in relationships by reference to historical landslide dams associated with earthquakes and heavy rain (61 disasters; 168 cases [44]). This theory quantifies the impacts of various factors on the several categories of dams and presents the findings in a graphical format. An example is shown in Figure 4. Correlations were apparent between the vertical axis (the dam length) and the horizontal axis (the dam height), and other factors are plotted on their individual axes. The results suggest that the various topographic conditions of mountain collapses (e.g., landslides and deep-seated collapses) influenced the deposition shapes of landslide dams.

Next, to simply predict the deposition heights and lengths of landslide dams based on only limited information (the “factors”), historical data on landslide dam formation were subjected to statistical analysis [44,45]. Specifically, using variables such as the topographic conditions shown in Figure 5, an optimal predictive formula was derived via repeated reduction and multiple regression analyses. Additionally, Figure 4 is expressed in dimensionless terms [56]. The conceptual diagrams that reveal the mountain collapse shapes associated with landslide dam formation (river channel blockages) are shown in Figure 6, and those that show the deposition shapes of landslide dams are presented in Figure 7 [45].

Additionally, the factors affecting the shape of the natural dam, as shown in Figure 5, are referenced. The formulae that predict the deposition height and length of a landslide dam, statistically derived via analysis of historical landslide dams [44], are as follows [45]:

where H_earthquake is the deposition height (m) of a landslide dam caused by an earthquake; H_{heavy rain} is the deposition height (m) of a landslide dam caused by heavy rain; L_earthquake is the deposition length (m) of a landslide dam caused by an earthquake; L_{heavy rain} is the deposition length (m) of a landslide dam caused by heavy rain; h₁ is the relative height from the top of the collapsed mountain to the riverbed (m); l₁ is the horizontal distance from the front soil movement point to the riverbed (m); w₁ is the damming width (the river channel width; m); A₁ is the area of collapse (/100,000 m²); B₁ is the width of the collapse (m); and α₁, β₁, γ₁, δ₁, α₂, β₂, γ₂, and ε₂ are the coefficients for each factor (Table 1) (Figure 6 and Figure 7).

To verify the reproducibility of our predictive formula for heavy rain, we used new data on landslide dams that formed in the Kii Peninsula of Japan during Typhoon 12 in 2011 (Talas) [50]. Figure 8 compares the predicted deposition heights and lengths to the actual values [50]. Figure 8 also shows the predictions for the historical landslide dams caused by heavy rain [44] that were used in the original analysis [45]; the predicted and actual values were well-correlated (correlation coefficients r: 0.75 and 0.91, respectively; generally, a correlation coefficient of 0.7 or above is considered to indicate a good correlation).

Next, we used Equations (1) and (2) to predict the deposition shapes of landslide dams caused by earthquake-induced embankment collapses that had not been previously analyzed, such as those in the Nigawa Yurino area (Table 1). Terrain and disaster data were retrieved from the literature [27,28,29] and the Nigawa Yurino area landslide museum (Table 2). Constructed embankments we viewed as components of the mountainsides. Constructed embankments were viewed as components of the mountainsides. The reproducibilities of the predictive formulae under such conditions (Table 2) are shown in Figure 9; the predicted and observed deposition lengths and heights in the Nigawa Yurino area were closely correlated (◆); the model was valid.

The historical topographical changes associated with embankment construction in the Nigawa Yurino area are shown in Figure 10 [57,58]. Historically, embankments spanned streams; they were thus large, valley-filling structures [19] that are generally considered to be at high risk of sliding collapse.

To consider how differences in the duration of water accumulation after river channel blockages and the time required for recovery and reconstruction [37] affected the disaster risk, we examined the risk of breach at the time of landslide dam formation. We evaluated many large landslide dams that developed in 52 locations (Figure 11; [58,59,60,61]) in the Imogawa River Basin (area: 38.4 km²) during the Chuetsu Earthquake, 2004. The water accumulations immediately after the earthquake (October 2004) are shown in Figure 12 (left), and those about 2 months later (December 2004) are shown in Figure 12. More than half of all landslide dams remained at full capacity 2 months after blockage. The rainfall in the basin during these 2 months, recorded by nearby AMeDAS stations [62], is shown in Figure 13.

As shown in Figure 13, the maximum hourly rainfall (upper bar graph, right legend: R mm/h) during this period was about 10 mm/h, and the cumulative rainfall over the 2 months post-disaster (lower line graph, left legend: ΣR mm) was approximately 120 mm. These results suggest that, in general, the risk of breach attributable to rainfall was low during this period. However, this varied by the basin area and the rainfall conditions. Thus, appropriate ameliorative measures are required within no longer than 2 months if landslide dams form under similar conditions (e.g., basin conditions) in the future. Special care is required in areas in which it is difficult to implement safety measures within the defined period (e.g., areas that large vehicles cannot easily access).

Figure 14 presents an overview of the disaster risk associated with river channel blockages triggered by the sliding collapse of embankments adjacent to rivers based on previous studies/reports. The following concerns are relevant: (1) direct damage to nearby residences when an earthquake triggers an embankment collapse; (2) landslide dam formation; (3) the inundation of upstream areas by accumulated water; and (4) flood and debris flow inundating downstream areas after landslide dams breach. However, these insights are based on past data; further exploration of, for example, topo-graphical conditions is desirable in the future.

The risks associated with flood outflow when landslide dams breach (point 4 above) are shown in Figure 14, and are discussed below.

Most breaches of landslide dams formed by river channel blockages can be classified into three types: erosion attributable to dam overtopping, large-scale collapse after water infiltration, and progressive collapse caused by piping [43]. In the past, breaches caused by overtopping erosion were the most common [42]. However, many breaches were not well-studied; further clarification was needed. Therefore, we performed several breach experiments using the terrains of real streams (for example, Figure 15) [54,55] in an effort to understand how landslide dams breached and the flood outflows from such dams. We verified the reproducibility of our experimental results using analytical models that predict such breaches [55]. However, we used only models of steep landslide dams (dam slope gradient θ_d: 1/2); gently sloping landslide dams (dam slope gradient θ_d: 1/3) require more attention (Figure 2 and Figure 3). We plan to further explore the overtopping breaches of landslide dams with gentle slopes.

The equations of water motions (flows) in the x-direction (downstream) and y-direction (transversely) in a two-dimensional field (the “shallow water flow equations”); the continuity equations for flow and sediment transport; and the continuity equation for the riverbed are as follows:

where x is the downstream distance, y the transverse distance, u the average flow velocity in the x direction, v the average flow velocity in the y direction, t time, ρ the density of the flow layer, p pressure, ε the eddy viscosity coefficient, τ_x and τ_y the bed shear forces in the x and y directions, h the flow depth, g gravitational acceleration, z the riverbed height, C the sediment concentration in the flow layer, C_* the deposition concentration at the riverbed, and i_b the rate of sediment erosion or deposition. The bed shear resistance was calculated using the Manning formula; we employed the hydrostatic approximation of pressure. The bed erosion rate was calculated using a traditional model that incorporates sediment concentration [63].

In contrast, to consider the widening of the water path due to lateral erosion, we used a model of the bank erosion velocity [64] that considers how secondary flows and local riverbed gradients affect the erosion rate. When the riverbed elevation differs between adjacent computational points, the effect of bank erosion is incorporated into the bed erosion rate i_bj of the higher riverbed elevation (computational point j) (Figure 16).

where Δz is the riverbed elevation difference between adjacent computational points and Δy the distance between these points. Note that erosion velocity is considered in four directions: upstream, downstream, left, and right.

Ashida et al. [65] found that when the dimensionless shear force acting on the bank was τ_*s, the bank erosion velocity was as follows:

Assuming that the friction velocity acting on the bank is proportional to the flow velocity near the bottom, u_s, the dimensionless shear stress τ_s acting on the bank can be expressed as follows:

where τ_c is the dimensionless, critical shear force; ρ the unit weight of sediment near the bottom; θ the internal friction angle; d the particle diameter; r the flow radius of curvature; U_L the flow velocity near the bottom; U the principal flow velocity; V the secondary flow velocity; and N and a₃ coefficients. Shimizu and Itakura [66] and Takahashi et al. [67] determined the radius of curvature r and other parameters.

To verify the reproducibility of this analytical model for an overtopping breach (for example, Figure 17), the predicted and observed (actual) values of the dam flood outflows at the time of breach are shown in Figure 18, which includes an analysis of both gentle (dam slope gradient θ_d: 1/3) and steeper (dam slope gradient θ_d: 1/2) [45] slopes. As shown in Figure 18, although some errors in flood peak timing are apparent, the results are generally robust.

Finally, using this analytical model [55], we predicted the flood outflow volume (the peak flow rate) from a landslide dam (with accumulated water) and the dam breach that triggered a river channel blockage in the Nigawa Yurino area. The prediction employed the known river channel topography and blockage shape [27,29] and topographical information [57,58] (Table 3). The analytical model predictions [55] made under the conditions of Table 3 were plotted against the international, empirical dam factor relationships (on a log-log scale) [42] and are shown in Figure 19. The peak outflow volume after the overtopping breach of a landslide dam that developed after an embankment sliding collapse in the Nigawa Yurino area was smaller than those of landslide dams caused by general mountain collapses but larger than those of glacial lake outburst floods [68,69]. This is partly explained by the fact that landslide dams that develop after embankment sliding collapses exhibit gentler deposition slopes (the deposition height is low compared to the deposition length) than general mountain collapses, given the mechanical properties of the embankment soil. This affects the dam outflow. The different behavior of glacial lake outbursts may be attributable to variations in the breach (melting) times, given that the erosion and thawing mechanisms are distinct, affecting the flood outflow. Further analyses of these and other cases are warranted.