Identifying post-earthquake debris flow hazard using Massflow

Catastrophic debris flows are common after large earthquakes and pose a significant risk for recovering communities. The depositional volume of these large debris flows is often much greater than the initiation volume, suggesting that bulking of the flow plays an important role in determining their volume, speed, and runout distance. Observations from recent earthquakes have driven progress in understanding the relationship between triggering rainfall events and the timing of postearthquake debris flows. However, we lack an adequate mechanism for quantifying bulking and applying it within a hazard context. Here we apply a 2D dynamic debris flow model (Massflow) that incorporates a process-based expression of basal entrainment to understand how debris flow bulking may occur within postearthquake catchments and develop hazard maps. Focussing on catchments in the epicentral area of the 2008 Mw 7.9 Wenchuan Earthquake, we first parameterised the model based on a large debris flow that occurred within the Hongchun catchment, before applying the calibrated model to adjoining catchments. A model sensitivity analysis identified three main controls on debris flow bulking; the saturation level of entrainable material along the flow pathway, and the size and position of initial mass failures. The model demonstrates that the difference between small and very large debris flows occur across a narrow range of pore-water ratios (λ). Below λ=0.65 flows falter at the base of hillslopes and come to rest in the valley bottom, above λ=0.70 they build sufficient mass and momentum to sustain channelised flow and transport large volumes of material beyond the valley confines. Finally, we applied the model across different catchments to develop hazard maps that demonstrate the utility of Massflow in post-earthquake planning within the Wenchuan epicentral region.


Introduction
Debris flows are water-saturated masses of soil and rock that rush down hillsides, funnel into stream channels, and evacuate large volumes of material from valley catchments (Iverson and George, 2014). They are a destructive natural hazard that threatens life and infrastructure in steep mountainous regions (Petley, 2012). Debris flow hazards depend upon a suite of processes governing the triggering of the flow and the entrainment of sediment and water during runout. Triggering processes affect the location, timing, and to a certain extent, the final volume and runout distance of a debris flow. The triggering mechanism of debris flows has received a lot of attention in the literature, and a number of possible mechanisms have been identified, such as shallow landsliding , Hortonian overland flow (Domènech et al., 2019), and entrainment of material during flooding (Pierson, 1982). After initiation, the characteristics of the flow pathway control the run out dynamics, (Braun et al., 2018;Cuomo et al., 2014;Fan et al., 2017). Of these characteristics, the entrainment rate of basal sediment is particularly important, as it affects the total volume of the flow, its runout length, and velocity (Berger et al., 2011;Frank et al., 2015;Pirulli and Pastor, 2012). Physical experiments have demonstrated that the rate of basal entrainment is non-linear and depends on the material properties of the bed sediment and their hydrology (Iverson et al., 2010). Measuring the hydrology and sediment properties along debris flow paths in the field is challenging (e.g. McArdell et al., 2007), and it is currently impossible to measure these properties in real time. This makes accounting for the non-linearity of basal entrainment a significant challenge, yet one that is necessary to develop effective debris-flow hazard models.
Estimating the potential extent of the largest debris flows is an important hazard challenge during the response and recovery phases after a large earthquake. Debris flow fans are often the only flat land within mountainous regions, but resettlement in these areas can increase exposure to debris flow hazards. Typically, empirical-statistical methods are employed to calculate hazard as they provide information on the timing of potential debris flows, can be readily implemented using limited data, and have potential for use in early warning (Huang et al., 2015). However, these methods do not provide information about the final volume of debris flow deposits. The statistical distribution of final landslide volumes is heavy-tailed, hence the amount of sediment that is deposited in these fans also follows a heavy-tailed distribution (Benda and Dunne, 1997). When a rainfall event triggers a debris flow, the volume of the resulting debris flow is not a simple function of the intensity or duration of the rainfall trigger. Numerical modelling (both process-and rule-based) has provided an avenue to understand the volume of debris flows (Frank et al., 2015). The runout model LAHAR-Z demonstrates how simple rules defining the relationship between debris flow volumes, widths, depths and topographic parameters can provide a method for estimating the spatial extent of debris flow hazards (Iverson, 2000); however, it requires prior understanding of the final depositional volume. Currently, few dynamic run out models include processbased descriptions of entrainment, but rather employ calibrated bulking factors that increase the flow volume proportionally to either topographic variables (slope angle, morphology, distance travelled) (Cascini et al., 2011;Fannin and Wise, 2001), or velocity and depth of the flowing mass (Chen et al., 2006;Cuomo et al., 2016;McDougall and Hungr, 2005). These approaches disregard changing local conditions that can drastically alter the final volume and run out extent.
Bulking of debris flows occurs non-linearly through a positive feedback where increased bed pore pressures cause scour of the bed and reduce basal friction (Iverson et al., 2010). The conditions under which this positive feedback occurs are poorly understood in a field context and are a significant challenge for hazard modelling. At the heart of the problem lies the heavy parameterisation required to adequately model the processes driving entrainment (Ouyang et al., 2015a). The hydrology of sediment along debris flow paths is crucial to the development of large debris flows (McCoy et al., 2012), as positive pore pressures can build in wet bed sediments, reducing intergranular friction and accelerating scour of the bed material Cuomo et al., 2016). The evolution of pore-water pressures during the propagation of debris flows largely determines the run-out distances and growth of the propagating masses . Within loosely compacted bed sediments, such as co-seismic landslide deposits, the rapid undrained loading of an overriding debris flow can lead to the sudden compaction of bed material, causing a surge in pore-water pressures (Pirulli and Pastor, 2012;Sassa, 1985;Sassa and Wang, 2007). Where there is sufficient water content, this surge can induce liquefaction of the bed material, accelerating entrainment, and causing a positive feedback that increases flow speed, mass, and momentum (Iverson et al., 2010;Reid et al., 2011). As such, relatively small variations in the initial conditions of the bed material can profoundly influence landslide behaviour . This paper investigates whether a process-based debris flow model can be used as a tool for estimating the likely extent of debris flow hazard to aid postearthquake recovery planning. We examine the relative importance of three factors that contribute to debris flow building: the size of the initial mass failure, the position of the initial mass failure, and the saturation of the entrainable material across the catchment. To assess the sensitivity of debris flow evolution to each of these controls, we explore the requisites necessary to initiate and sustain the largest events using a 2D dynamic model of debris flow evolution (Massflow) that incorporates a process-based expression of basal entrainment (Ouyang et al., 2015a). We then apply this to a post-earthquake context by assessing three catchments affected by the 2008 M w 7.9 Wenchuan earthquake.

Methods
Massflow is a depth integrated mass and momentum conservation model that simulates flowing water-sediment mixture over an erodible bed layer (Ouyang et al., 2015a). It incorporates a formulation of basal entrainment that is a function of flow and bed material properties. As an initial calibration, and to ensure model consistency, we compared our model outputs against a prior application of Massflow at Hongchun Gully (Wenchuan, China) by Ouyang et al. (2015a), and assessed the impact of a modification made to limit the depth of entrainment to reflect the limited co-seismic landslide deposits available within the catchment. We then used this process-based expression of bed entrainment to explore the potential debris flow hazard by performing a sensitivity analysis to examine the relative importance of initiation location, volume, and bed hydrology in controlling the size and final run out extent of debris flows. Finally, we applied our model to predict debris flow potential in nearby catchments triggered during the same intensive rainfall event and produced hazard maps.

Study areas
We focus our analysis on the epicentral area of the Mw 7.9 2008 Wenchuan earthquake in the Longmen Shan, China. The Longmen Shan are a transitional mountain belt between the Sichuan Basin and the Western Sichuan Plateau, characterized by rugged mountains interspersed by deeply incised valleys ranging from 860 to 3950 m elevation. The Longmen Shan Fault Zone, located on the southern border of our study area, trends northeast southwest and generates large earthquakes (Burchfiel et al., 1995). Granitic rocks, Sinian pyroclastic rock, Carboniferous limestone, and Triassic sandstone underlie the area, with loose Quaternary deposits distributed along terraces and alluvial fans (Tang et al., 2011).
The 2008 Wenchuan Earthquake mobilised an estimated 3 km 3 of material across the Longmenshan (Li et al., 2014). The co-seismic landslide deposits provide sources for rainfall induced debris flows. For example, a single intensive rainfall event spanning the 13th to 14th of August 2010 triggered > 20 debris flows in catchments along the Min Jiang (Tang et al., 2011); among these were Hongchun, Bayi, and Yinxingping (Fig. 1). These fifth order or larger catchments (> 5 km 2 ) are in close proximity. All catchments have steep (30°-60°) valley walls and high relief (> 1 km), and all are underlain with the same Precambrian granitic bedrock. Axial channels are steep (~20°) and empty into alluvial fans at the valley mouth (Table 1). Each catchment also contains large volumes of co-seismic landslide deposits (> 3 × 10 6 m 3 ), yet only Hongchun and Bayi produced debris flows that transported large volumes of material onto the alluvial fan at the valley outflow.

Hongchun
Hongchun gully is located along the Min Jiang, on the opposing bank to Yingxiu town, which was devastated during the Wenchuan earthquake in 2008. In 2010, a debris flow from the Hongchun gully partially dammed the Min Jiang, diverting its flow and flooding the town. The Hongchun catchment covers an area of 5.35 km 2 ; the main channel length is 3.55 km with an average gradient of 19.5° (Table 1). The August 14th 2010 rainfall event triggered three simultaneous debris flows in each of the branch gullies, contributing initial source volumes of 11.2 × 10 4 m 3 , 3.9 × 10 4 m 3 , and 3.2 × 10 4 m 3 , which combined in the main channel entraining material along the valley bottom. The final volume of material deposited on the alluvial fan at the valley mouth was~80.0 × 10 4 m 3 , with~40.0 × 10 4 m 3 being carried into the Minjiang River (Xu et al., 2012).

Bayi
Bayi gully is situated near Longxi town on the banks of the Longxi River. The catchment covers an area of 8.3 km 2 and its western edge borders the Hongchun valley catchment; the main channel length is 4.23 km with an average gradient of 21.3° (Table 1) (Chang et al., 2017). The final extent of the large debris flow triggered by the rainfall event of August 14th is clearly discernible in post-event satellite imagery, with measurements of deposition depths up to 7 m.

Yinxingping
Yinxingping is located 14 km north of Hongchun gully along the western bank of the Min Jiang. It has a catchment area of 7.06 km 2 ; the main channel length is 4.41 km with an average gradient of 24.1° (  Table 1). Despite large co-seismic landslide deposits within the valley basin (10.3 × 10 6 m 3 ), and numerous debris flows initiated within the valley confines, the final run out extent was small by comparison to the other catchments, with an estimated 3.0 × 10 4 m 3 deposited at the valley mouth (Tang et al., 2011).

Massflow model description
We use a modified version of the Massflow model originally developed by Ouyang et al. (2015aOuyang et al. ( , 2015b. It employs a second order solution of the depth integrated mass and momentum balance equations (shallow water equations) to obtain flow heights and velocities across a Cartesian grid inclined parallel to the pervading slope, modified to include the momentum exchange across the boundary layer due to entrainment (Tai and Kuo, 2008). It is a single-phase continuum model that routes a mass of mixed solid and fluid particles over an erodible base layer while incorporating basal entrainment at a rate (E) defined as where z b is the elevation of the debris flow base, ρ _ is the depth averaged density of the flowing mass, u _ and v _ are flow velocities in the x and y directions, τ 1b is the basal traction of the flowing mass, and τ 2s is the resistive shear stress from the erodible base layer. The basal traction of the flowing mass is expressed as a coupled Voellmy-Coulomb friction model to capture increased traction at higher velocities without requiring an internal friction angle smaller than the expected angle of repose along the deposition surface (Ouyang et al., 2015a), such that ; φ voellmy is the internal friction angle of the flowing mass, C z is the Chezy coefficient, ρ w is the density of water, h is the flow height, g z is the component of gravity acting normal to the inclined slope, and δ is the basal friction angle. The resistive shear stress of the base material is assumed to follow the failure criterion: where c and ϕ 2 are cohesion and friction angle of the erodible material, and λ is the pore-water ratio that indicates the degree of saturation of the bed material. For a full derivation of the entrainment function, see Ouyang et al. (2015a).
To solve the shallow water equations we use a MacCormack total variation diminishing (TVD) finite difference scheme (Liang et al., 2007;Ouyang et al., 2013;Ouyang et al., 2015b). To simulate the  limited material available across a bedrock landscape, we impose a maximum entrainment depth homogeneously across the landscape. Any deposited material that subsequently raises the elevation above the limit is afterwards made available for entrainment.

Model calibration and consistency
As an initial verification of model consistency and an assessment of the impact of our modification, we replicated the August 14th 2010 Hongchun gully debris flow using the initial failure conditions and parameterisation described in Ouyang et al. (2015aOuyang et al. ( , 2015b (Table 2).
We modified Massflow to specify a maximum depth of entrainment to limit unrealistic debris flow bulking. We set the limit of entrainment to 6 m consistent with the final deposition depth at the Hongchun alluvial fan, adjusting the parameter incrementally until our results corresponded with those of Ouyang et al. (2015aOuyang et al. ( , 2015b. Having calibrated our model formulation, we performed the simple sensitivity analysis described in Ouyang et al. (2015aOuyang et al. ( , 2015b to evaluate the effect of imposing a limit on entrainment and the coarser resolution of our simulations. We applied Massflow across a 20 × 20 m Cartesian grid overlain by a globally available digital elevation model of the same resolution, as higher resolution data was not available. The application of this model to a globally available DEM, demonstrates the potential to transfer the model to other settings. Initially, we simulated the simultaneous failure of landslide deposits at the head of three branch gullies (Ganxipu gully, Dashui gully and Xindianzi gully) with volumes of 11.2, 3.9, and 3.2 × 10 4 m 3 respectively. We ran the simulations until such time that all material within the valley confines were static (i.e., no longer being entrained, transported, or deposited).

Sensitivity analysis
To explore the sensitivity of debris flow bulking to small variations in environmental conditions we undertook a sensitivity analysis to determine the main control limiting the final volume of extreme events. We systematically varied each of the parameters expressed within the entrainment function that are explicitly subject to change between events. These are limited to the thickness of the initial failure volume (h), the position of the initial failure volume, and the pore water ratio of the bed material (λ) -though the frictional properties of the bed material may vary implicitly as a result of varying water content. We concentrated our analysis on the Hongchun gully debris flow, varying the size of the initial failure volume, the position of the initial failure volume, and λ for a set of expected values to produce scenarios that encompass the range of possible extreme events.
We assume that three simultaneous failures triggered debris flows at Hongchun, each at the head of branch gullies. To assess the relative importance of initiation volume we scaled the observed source volumes from a total of~20 × 10 4 m 3 to 5, 10, 30, and 40 × 10 4 m 3 whilst maintaining their original positions. To assess the impact of initiation position on the final volume and extent of debris flows, we released single source volumes of 5, 10, 20, 30 and 40 × 10 4 m 3 from the head of the Ganxipu branch gully at positions ranging from 150 to 650 m of vertical relief on the steep slopes measured as elevation above the valley channel. We modelled each combination of initiation volume and location across the expected range of pore water ratios (λ: 0.5-0.8), enabling a direct comparison of all three variables to assess the dominant control on determining the magnitude of debris flow bulking.

Translation of calibrated model
We attempted to reproduce two debris flow events in Bayi and Yinxingping by directly applying the configuration calibrated for Hongchun, varying only the topography and initiation source volumes. We identified debris-flow source areas from post-event satellite imagery, and used an empirical area-volume scaling relation derived from a landslide inventory of the study area (Tang et al., 2011) to estimate initiation source volumes used as model inputs: where D is average landslide depth, S L is landslide area.

Model results
The main difference between the original formulation of the Massflow model in Ouyang et al. (2015aOuyang et al. ( , 2015b and our version is the implementation of a limit on the depth of scour. This modification affects the spatial distribution of modelled erosion and deposition and the final volume of material deposited beyond the valley confines onto the alluvial fan (Table 3). However, comparable sensitivity analyses conducted by Ouyang et al. (2015aOuyang et al. ( , 2015b) that altered the pore water ratio (λ) and cohesion (c) values result in larger deposited volumes than our model (Table 3). This is most evident where the pore water ratio is raised to 0.75, where modelled deposition volumes differ by 2.4 × 10 6 m 3 . The difference relates to limited scour depth in our modified model. In the original formulation and average of~12 m entrainment occurs along the entire length of the valley bottom, an observation that cannot be confirmed by field observation. Also, the total volume of the deposited material is larger than the total estimates of coseismic landslide deposits within the catchment (Tang et al., 2011;Xu et al., 2012), suggesting that the simulated debris flow would need to entrain non-earthquake derived sediment during bulking. Hence, setting a maximum depth of entrainment enables an assessment of the likely volume attained by the largest events limited by the material available across the landscape, whilst still capturing the variability of debris flow bulking between parameter sets.
The sensitivity analysis demonstrates that the final volume of material deposited on the alluvial fan grows by an order of magnitude across the expected range of each parameter (Figs. 2-4). Changes to the vertical height above the valley channel of initial mass failures, while maintaining a pore-water ratio of 0.7, show that landslides that create debris flows that exit the catchment either have large initial failure volumes, or initiate close to ridges. Of the two parameters, the size of the initial failure volume acts as the stronger control on final debris flow size (Fig. 2). We think that this reflects the balance between the velocity of the flow that enters a channel system (the relief factor) and the magnitude of the basal shear applied by the flow (the volume Table 2 Model parameters for Hongchun Gully (Ouyang et al., 2015a(Ouyang et al., , 2015b.  factor). At elevations above the channel between 0 and 300 m, the initial mass failure attains low velocities and is unable to entrain material to the imposed limit. Increasing the size of the initial failure volume at these low elevations incrementally increases the basal shear force, and so the rate of entrainment along the flow pathway, but not to the point where the imposed limit is reached. Flows that initiate > 300 m above the channel network attain velocities high enough to entrain channel sediment. Hence flows with smaller initial failure volumes that initiate close to ridges have the potential to create valley exiting debris flows through bulking within the channel network.
Comparing the size of initial failure volume (Fig. 3) and the relief (Fig. 4) factors with the variability in the pore water ratio in the bed demonstrates the sensitivity of the model to bed pore water. The rate of change in the final deposit volume with elevation slows above initiation volumes of 15 × 10 4 m 3 when maintaining a constant initiation elevation (515 m of relief) and varying the pore water ratio of the bed material ( Fig. 3b). At the highest pore water ratio (0.75), the proportional increase in debris flow magnitudes for initiation volumes above 20 × 10 4 m 3 are much reduced, as the limit of entrainment is reached along the entire flow pathway and extending the area scoured by the flowing mass is hampered by the confines of the valley walls. Below a pore-water ratio of~0.65 there is little to no material deposited at the valley mouth, though once exceeded final deposition volumes increases rapidly (Fig. 3c). The growth of final deposition on the alluvial fan relative to pore-water ratio shows little distinction between the size of initial failure volumes, with initial failure volumes of 30 and 40 × 10 4 m 3 resulting in debris flows separated by little more than the initial difference (Fig. 3c). This again suggests that once the initial failure volume is sufficient to entrain material along the flow path to the depth limitation, additional growth is impeded. By contrast, the final volume of deposited material relative to initiation volume are vastly different between pore-water ratio scenarios, indicating that the water content of the bed material is the dominant constraint upon the final volume of material deposited on the alluvial fan (Fig. 3).
Maintaining a constant initiation volume of 20 × 10 4 m 3 and varying both the elevation of initiation and pore-water ratios, we see the same properties and limitations as in the previous analysis. Increasing initiation elevation linearly increases the final magnitude of  the debris flow, but with considerable variation between pore-water ratio scenarios (Fig. 4b). Where pore-water ratios are below 0.65, there is no deposition of material on the alluvial fan at the outflow of the valley irrespective of the relief of the initial mass failure. Unlike initiation volumes, the magnitude of final deposits do grow with higher elevation scenarios as the extended flow pathway provides additional material for entrainment, though the increase is comparatively small (Fig. 4c).
The results of our sensitivity analysis demonstrate the dominance of pore water conditions in controlling the volume and run out extent of the largest debris flows. The location and initiation source volumes provide less of a control on whether a debris flow will be large enough to exit a fifth order catchment. These results provide a framework for assessing the hazard posed by the largest debris flows in a particular catchment. An example of where this may be particularly important is during the recovery phase after a large earthquake, when resettlement or recovery activities may occur on alluvial fans. Hence, a first order estimate of large debris flow risk can be obtained even with poorly constrained initiation locations and volumes. However, the pore water content of the material entrained in the flow ultimately controls the volume of the final flow. Moving beyond estimates of potential hazard into real-time warnings of large debris flows requires data on bed pore pressures that is currently impossible to collect in the field.

Hazard applications
To test the utility of our modelling approach for hazard analysis, we applied the calibrated model in two different ways. Firstly, we calculated the extent of large debris flows in our two other catchments associated with August 2010 rainfall event. Secondly, we created maps that show the extent of plausible debris flow sizes for the same catchments, and use this to create a tool for better spatial planning of potential hazards.

Hindcasting debris flows associated with the August 2010 Yinxiu rainfall event
We applied our calibrated model at both Bayi and Yinxingping to hindcast the large debris flows observed during the August 2010 rainfall event. We identified the location of initial failure volumes by analysing pre-and post-event satellite imagery and approximated the size of the initial failure volumes using an area-volume scaling relation (Eq. (4)). We then assumed these initial source volumes failed simultaneously and simulated the evolution of the resultant debris flows. We assumed that the bed pore-water ratio was the same as that calibrated at Hongchun.
The Bayi catchment produced a large debris flow during the 2010 event. We modelled the simultaneous triggering of seven landslides at the head of branch gullies within Bayi. The run out extent of the resultant debris flow is in close agreement the mapped deposit extent (Fig. 5). We calculated the final deposit volume to total 125.9 × 10 4 m 3 , which is consistent with field investigations that report deposition depths of up to 7 m, and a total of 116.5 × 10 4 m 3 deposited on the alluvial fan (Ma and Li, 2017).
Yinxingping catchment provides a contrast, as there was only a small debris flow event initiated during 2010. We modelled the simultaneous triggering of four mass failures at the head of branch gullies discernible from post event satellite imagery. Yinxingping produced small debris flows that combined in the valley channel, but did not coalesce into a large debris flow, with observations reporting deposition depths of up to 3 m on the alluvial fan, with a total volume of 3.0 × 10 4 m 3 (Tang et al., 2011). However, our model formulation drastically overestimates the final run out extent and volume of the resultant debris flow, predicting depths of up to 12 m and a final deposition volume of 62.3 × 10 4 m 3 (Fig. 6).
The contrast between the results of the Bayi and Yinxingping hindcasts reflects the challenges of model parameterisation (particularly the bed pore-water ratio) even where data is from an adjacent catchment for the same storm event. The contrasting results also highlight the limitation of this model as a tool to forecast the magnitude of individual flows associated with specific rainfall events. That the model failed to replicate the debris flow observed at Yinxingping using the parameterisation of Hongchun is unsurprising given that one or more of the controlling characteristics differed between the two catchments. Possible differences in topography, properties of the bed material, the availability of entrainable material, and the water content of the bed material (pore-water ratio) may be the cause of these discrepancies. Our models of topography (DEM) at Yinxingping are the same as in Hongchun and Bayi, so are unlikely to have accounted for the differences in deposit volumes. Our estimates of the volume and location of co-seismic deposits, which provided the source material for debris flow initiation in both catchments are also well constrained using remotely sensed imagery. While the mapping of the size and location of landslides is accurate, it is difficult to account for the material Fig. 4. A) Final deposition volumes (colour scale on RHS) relative to both the pore-water ratio of the bed material (λ), and the elevation of the initial failure volume (v). B) Final deposition volumes plotted against the elevation of the initial failure volume for pore-water ratio scenarios (0.6-0.75). C) Final deposition volumes plotted against the pore-water ratio for initial elevation scenarios (145-645 m).
properties of each landslide and how these have evolved since the earthquake (Domènech et al., 2019). From the satellite imagery (resolution of~1 m) the style and distribution of large grains is similar for the largest landslide in each catchment. The pore water conditions at the bed are another likely candidate for these differences. TRMM data during the 48 h spanning the events shows that Hongchun received 37% more cumulative rainfall than Yinxingping; 148 mm compared to 108 mm. In addition, the higher relief, and steeper terrain and valley channel gradient at Yinxingping may have facilitated faster draining of the catchment, reducing the saturation level of the confining walls and valley bottom relative to Hongchun. Therefore, the inability of the model to accurately forecast the size and extent of the debris flow at Yinxingping could be the result of the disparity between the pore-water ratios at the two locations. Although we did not undertake a formal calibration at Yinxingping, whilst conducting the hazard assessment (presented in the following section), we found that reducing the porewater ratio to 0.65 simulates debris flows that more closely match those observed, which supports our assertion that the entrainable material was less saturated at Yinxingping than either Hongchun or Bayi.
By contrast, the successful implementation of the model at Bayi may reflect the similar drainage properties and rainfall characteristics of the adjacent Hongchun valley. Local rain gauges describe similar patterns of precipitation at the base of Hongchun and Bayi, with total cumulative rainfalls across debris flow events of 143 and 124 mm respectively (Zhou and Tang, 2014). In addition, the catchment relief, and channel slope are in closer agreement (Chang et al., 2017;Tang et al., 2011), suggesting that these two catchments drain within the same timeframes, making it likely that the pore-water ratio of the bed material was consistent between these two sites.
These examples exhibit the challenges associated with hindcasting specific flows. Even in our relatively well-constrained locations, small changes in bed pore-water ratios may affect significant changes in debris flow volumes. As such, this modelling approach is not necessarily well adapted to predicting the size of a specific debris flow event.  However, it may be useful as a tool for understanding the range of possible debris flow event sizes to form the basis of preliminary analyses for categorising the most hazardous catchments.

Developing hazard maps
The sensitivity analysis performed on Hongchun gully demonstrated that the model is most sensitive to bed pore water ratio. While this is difficult to constrain in real time, we have a robust understanding of the ranges of bed-pore conditions that are possible in nature (Ouyang et al., 2015a). Therefore, we argue that this model is useful for planning; in particular assessing the potential for and possible extent of a large, catchment exiting debris flows. Hence, it is possible to define a range of debris flow scenarios if (1) we have an understanding of the approximate location and volume of initiation, and (2) if we can calibrate the model material parameters using a previous event. Our method contrasts other landslide-debris flow hazard assessments, where the probability of landslide failure is based on a quantification of static factors (slope angle, lithology, soil thickness, permeability, flow accumulation area, etc.), combined with a calibrated run out model (Chang et al., 2017;Stancanelli et al., 2017;van Westen et al., 2006).
Our methodology is to create maps of hazard potential by stacking the result of model predictions where we systematically vary the pore-water ratio (Fig. 7). To implement this modelling approach for our 3 test catchments, we first performed the local calibration using the extent and volume of the previous Hongchun event. Then, we chose the distribution of initial failure volumes by identifying the largest landslides that are close to the channel network and assumed they all failed simultaneously. We did not evaluate the likelihood of any individual landslide failing, the probability of simultaneous failures, or the effect of reduced pore-water pressures on initiation, but rather assume large landslides will occur and set up conditions to simulate the largest debris flow events that may be expected across a range of λ (0.5-0.8 with steps of 0.05). The effect of this assumption is that our modelled extents approximate the worst-case scenarios that we might expect each catchment to produce with increasing levels of saturation. For all runs, we keep material parameterisations and initiation volumes constant and delineate the extent of debris flow inundation above 2 m in depth. Areas that inundate across the entire range of pore-water ratios are assigned the highest hazard classification (Magnitude 7), reducing one level in severity for each scenario that fails to inundate the area, to the minimum classification (Magnitude 0) for areas that never inundate > 2 m.
The potential utility of these maps is demonstrated for our three study catchments in China. At both Hongchun and Bayi the zonation of high hazard restricted to the steep slopes of the valley walls and along the channel immediately below branch gullies ( Fig. 7a and b), whereas areas of high hazard extend slightly further down the main channel at Yinxingping (Fig. 7c). The actual debris flows recorded in the Hongchun and Bayi events represent a magnitude 3-4 hazard, while the magnitude 4-5 scenario best fits the Yinxingping event in 2010. When seen in the context of multiple runs, the difference in bulking associated with only small changes in pore water ratio may explain the dramatic differences in the size of specific debris flows (as in Section 4.1). In addition, there are some interesting patterns when examining all three catchments together. Where pore-water ratios are below 0.65 (magnitude 5-7), landslides triggered on hill slopes come to rest in the valley bottom but do not exit the catchment. However, where pore-water ratios are above 0.7 (magnitude 1-4), even small source volumes that fail on hill slopes high enough to build sufficient momentum to propagate along the valley channel, have the potential to entrain large volumes of material and pose severe risk to nearby inhabitants. Rather than a continuum of increasing debris flow magnitude inundating the valley outlet with increasing pore-water ratios, our analysis indicates a threshold that demarcates the volume of material deposited on the alluvial fan. Below the threshold, initial mass failures build as they scour the steep slopes, but immediately deposit material when encountering the valley channel, unable to transport the load at the reduced gradient. Above the threshold, the initial mass failures build sufficiently on the steep slopes to continue their momentum along the shallower gradient of the valley bottom and carry the load to the valley outlet. Whilst increasing the size and location of the initial mass failure extends the run out of the initial landslide along the channel bottom, only increasing the pore-water ratio above the observed threshold allows the largest debris flows to develop, and evacuates the accumulated material out of the valley confines. The resultant debris flows are bimodal in character, categorised by those that do, and those that do not, propagate along the valley channel. Our analysis identifies a narrow range across λ as being the main determinant of whether a particular debris flow will traverse a channel network to become a large debris flow (0.65-0.7). However, this is a site-specific threshold and subject to the material parameterisation imposed by Hongchun. The set of conditions suitable for sustaining valley channel flow will vary between catchments as the material properties of the bed sediments, spatial distribution of water content, and valley topography differ.
Finally, our tool may be used to better inform mitigation practices. For example, in China, the engineering standard for the development of debris dams relies on an equation describing debris flow discharges as a linear function of precipitation (Xu et al., 2012). Massflow provides an open source tool for the physical modelling of debris flows. For this model to be effective, we have demonstrated that it requires a sitespecific analysis of material properties obtained by back-calculation from an existing debris flow event. However, once calibrated it is possible to simulate the release of source volumes at the head of branch gullies and estimate debris flow volumes and spatial susceptibilities. The challenge of the back-calculation method is understanding the physical meaning of the pore-water ratio term. Currently, work needs to be done to understand how these ratios relate to measurable catchment properties, such as such as precipitation, temperature, grain-size distribution, depth to bedrock, drainage area, relief, slope angles, and valley channel gradient. Despite this challenge, it is possible to assign hazard levels based on the recurrence interval of conditions required to sustain flowing masses along valley channels. The final model produces a set of realistic volumes and velocities that can better inform mitigation measures. In addition, the tool provides a useful boundary object for the development of adaptive management strategies for potential landslide hazard.

Summary and conclusions
Between August 13th and 14th 2010, an intense rainfall event triggered > 20 debris flows in tributary catchments along the Min Jiang in the Wenchuan area of China (Chang et al., 2017;Tang et al., 2011). The numerous landslides triggered at both Hongchun Gully and Bayi Gully coalesced in their respective valley channels to form large debris flows that deposited huge volumes of material on alluvial fans at the valley out flows. Yet, despite the similarity in catchment characteristics, the majority of landslides initiated within Yinxingping Gully came to rest in the valley channel, with negligible material evacuating the valley confines. Here we attempt to replicate the observed debris flows at each of these catchments using a dynamic 2D numerical model of debris flow evolution (Massflow) that incorporates a process-based description of basal entrainment. We explore the conditions necessary to initiate and sustain the largest debris flows and assess the relative importance of key physical factors in controlling rates of basal entrainment to identify the root cause of disparities in debris flow magnitudes between adjacent catchments. We find that the pore-water ratio of erodible materials primarily controls the final volume and extent of debris flow run out, and that the size and location of initial mass failures is of secondary importance. We also find that, within our study areas, model scenarios that simulate mass failures at the head of branch gullies across a range of pore-water ratios (0.5-0.8) display evidence of a threshold, demarcating large debris flows that evacuate material out of the valley confines from smaller debris flows that come to rest in the valley bottom. Above the pore-water threshold of 0.7, sufficient water is available for debris flows to gather enough mass and momentum to continue entraining material along the length of the valley bottom and deposit large volumes of material on the alluvial fan. Below the porewater threshold of 0.65, the resistive forces along the valley channel overcome the accumulated momentum of the initial mass failure bringing it to rest before the valley outflow, and so restricting the hazard to the main valley channel.
Identifying the main drivers of debris flow bulking and constraining the conditions necessary to develop and sustain the largest events may aid future hazard identification and inform adaptive management practices that focus on prohibiting the conditions required for channelised flows (such as soil moisture reduction), rather than hard engineering solutions.