Experimental and Numerical Mudflows Modeling for Runout and Deposition Height Assessmen

doi:10.21203/rs.3.rs-3338795/v1

Download PDF

Research Article

Experimental and Numerical Mudflows Modeling for Runout and Deposition Height Assessmen

https://doi.org/10.21203/rs.3.rs-3338795/v1

This work is licensed under a CC BY 4.0 License

You are reading this latest preprint version

Mudflows are natural phenomena starting from landslides and presenting high impact when they occur. They generate great catastrophes in their path because most of the time there is no indication prior to the failure that triggers them. Understanding how mud is transported is of great importance in infrastructure projects that coincide with hillside areas due to the high risk of occurrence of this phenomenon by cause of the high slopes, which can involve great risks and produce disasters that involve great costs. This work presents the evaluation of mudflows, from the implementation of a laboratory scale experiment in a consistometer with its calibration and validation from numerical models to estimate rheological parameters of the material. Tests were also carried out in an open channel in the laboratory, based on the data previously obtained considering the behavior of the material as a both Newtonian fluid and non-Newtonian fluid. The experiment considered a channel with dimensions of 3 m long, 0.5 m high and 0.7 m wide with slope control, and a mud composition of silty material with 60% moisture. The tests were conducted with slopes of 5%, 10%, 15% and 20%. The numerical models were carried out in ANSYS FLUENT software. In addition, the calibration data of the numerical model were used for a real case study, simulating the slip flow occurred in Yangbaodi, in the southeast of China, occurred on September 18, 2002. The results of the numerical models were compared with the experimental results and show that these have a great capacity to reproduce what is observed in the laboratory when the material is considered as a non-Newtonian fluid. The model reproduced in an appropriate way the movement of the flow at laboratory scale, and for the aforementioned case study, some differences in the final length of deposition were noticed, achieving interesting results that lead the use of the calibrated model towards the estimation of risks due to the mudflow occurrence.

ANSYS

experimental test

landslides

mudflows

numerical simulation

rheological behavior

runout

Mudflow is a kind of debris flow that has been classified as a non-Newtonian flow, with large velocities and long run-out distances, which is excited by heavy rain in mountainous valleys (Bokharaeian et al. 2022; Qingyun et al. 2022; Pastor et al. 2023). Mudflows are natural phenomena starting from landslides that can have the potential to evolve into flows and thereby gather more channel sediments and cause significant destruction along its path and in the downstream areas (Ortiz-Giraldo et al., 2023),presenting a high impact when they occur, generating great catastrophes in their path, because most of the time there is no previous indication of the failure that triggers them (Li et al. 2012). Generally, mudflows are characterized by rapid movement, sudden occurrence, fluid-like movements, and long travel distances from their exit. Some studies indicate that the occurrence of mudflow is triggered mainly by excess pore pressure or liquefaction processes in the material, generating sudden and high shear strength losses in the landslide zone (Hungr 1995; Hungr et al. 2001; Lacerda et al. 2004; Zhuang and Peng 2014). Thus, the loss in shear strength occurs in parallel with a total fluidization of the soil, generating a mudflow.

In hillside tropical areas, landslides, and slope stability, as well as post-failure consequences that trigger mudflows, are issues of great importance in urban planning. In Colombia, high impact mudflow phenomena have been generated, for example, in the municipalities of Salgar and Mocoa, (Vega and Hidalgo, 2023; Hidalgo and Vega 2021; Vargas-Cuervo et al. 2019) at western and south of the country respectively. In Brazil, a high number of cases have been registered (Cabral et al. 2023). A large number of infrastructure projects and urban growth coincide with these areas, which normally have high slopes, and combined with the occurrence of high intensity rainfall and subsequent infiltration, lead to a greater probability of occurrence of these phenomena (Hungr et al. 2001; Lee and Widjaja 2013), generating high risks for people and high costs to cover the disasters ocurred, including the relocation of communities, repairs and reconstruction of the affected physical structure and restoration of supply networks (Vega and Hidalgo 2016; McDougall 2017; Xia and Liang 2018).

In spite of all the efforts in stabilization processes and slope prevention works to avoid landslide flows, this phenomenon continues to be one of the problems that most reorient economic resources for the repair of damages caused, and depending on the magnitude, it can claim human lives (Lacerda et al. 2004; Zhuang and Peng 2014; Vega and Hidalgo 2016; Vega and Hidalgo 2017). This confirms the need to continue implementing and searching techniques to improve the understanding of the behavior of slopes, constituting an interesting challenge that would allow to evaluate the risk to which future infrastructure projects would be exposed when they are located in areas threatened by the study phenomenon, and to prevent or mitigate the damage of existing infrastructure.

For the understanding of the behavior of slopes and its flow dynamics, modeling processes have been widely used, providing an alternative way to reproduce the behavior of a mudflow under the possible physical conditions under which the flow would be subjected, and therefore to quantify the outflow distance and velocity of the material, facilitating risk assessment and management (Iverson 1997; Cascini et al. 2010). Understanding the behavior of mudflows and accurately assessing the associated risks are crucial for effective hazard management and mitigation strategies.

In recent years, advancements in experimental testing joined to numerical modeling techniques, particularly employing Computational Fluid Dynamics (CFD) and Smoothed Particle Hydrodynamics (SPH), have provided valuable tools for studying and predicting mudflow behavior. The application of numerical methods in the rheological simulation of flows has developed rapidly because of advances in computing resources and numerical technology (Qingyun et al. 2022). Some works have used depth-averaged mathematical models to simulate all stages of debris flow, from its initiation to its disposal, using the finite volume method to solve the equations governing the flow behavior (George and Iverson 2014; Liu et al. 2018b; Xia and Liang 2018). This also makes it possible to determine the impact on the existing physical infrastructure, presenting an approach that could be implemented in the quantification and estimation of landslide risks. Some results have been successful and they have made it possible to reproduce landslides as case studies, showing the potential that the application of the models could have (Li et al. 2019; Feng et al. 2019; Peruzzetto et al. 2020).

Some research aims to study the mudflow as a heterogeneous medium, using numerical modeling processes from the Discrete Element Methodology (DEM), separating the solid particles from the fluid medium and obtaining a behavior of the interaction of these particles (Li et al. 2012; Liu et al. 2018a; Li and Zhao 2018; Zheng et al. 2018; Trujillo-Vela et al. 2020). In general, the results were positive in terms of predicting the behavior of a landslide and associated flows and the possible effects on existing structures. However, when they were coupled with other modeling procedures such as the Material Point Method (MPM), the results were positive to predict the behavior on landslides and the possible effects on existing structures (Liu et al. 2018a) or the use of CFD computational fluid dynamics (Li and Zhao 2018; Zheng et al. 2018).As a complement to simulate the solid and liquid part of a landslide-derived flow, DEM becomes a very effective method in granular type flows.

CFD applications are often used to model various flow problems and phenomena that sometimes cannot be easily described from simple conceptual models (Versteeg and Malalasekera 2007). This numerical modeling alternative is a tool that through its solution can contribute to the simulation of a mudflow when considered as a continuous medium. The study carried out by Trewhela et al. (2014) have been able to simulate the behavior of mudflow in pipelines. Moreover, the studies conducted by George and Iverson (2014); Baggio et al. (2021); Peruzzetto et al. (2022); Abraham et al. (2022); and Wang et al. (2023) have modeled debris flows in natural and artificial open channels, analyzing variables such as flow velocity and volume fraction, with the aim to predict the runout and possible effects of their final deposition.

Despite attempts to simulate the behavior of mudflows, even in real case studies, there is still a great shortage of models with CFD computational fluid dynamics solutions, where aspects such as viscosity from the water content of the soil, control the way the flow moves (Lee and Widjaja 2013; Buiskikh 2015; Widjaja and Parahyangan 2015) and its interaction with the physical and geometrical conditions of the surface over which it flows. Rheological studies of mud at sufficiently high concentrations of solids have shown it to be a highly viscous non-Newtonian fluid, exhibiting a creep limit evidenced by an observed minimum depth necessary for a uniform layer of mud to flow (Huang and García 1998; Huang & Aode, 2009; Carotenuto et al. 2015; Xu and Huhe 2016; Kameda and Hirauchi 2018). Hence, there are rheological models that simulate the behavior of the material as a viscous flow once established within the model, which can be coupled to simulate the flow advance (Sawyer et al. 2012; Widjaja and Lee 2013; Han et al. 2019).

Therefore, by combining the insights gained from experimental tests and the computational capabilities of CFD, researchers and engineers can develop comprehensive rheological models that offer detailed insights into mudflow dynamics, deposition patterns, and potential extents of propagation. By calibrating numerical models with experimental data, the accuracy of the simulations can be enhanced, ensuring realistic representations of mudflow behavior. Some examples of this kind of approach can be found in Guo et al. (2023) and Bokharaeian et al. (2022). This approach allows for the examination of various scenarios by manipulating input parameters such as moisture content, slope steepness, and flow types. Consequently, it becomes possible to evaluate the influence of these factors on the initiation, flow dynamics, and deposition patterns of mudflows. By simulating mudflow events in specific areas prone to landslides, it is possible to estimate the potential impacts on infrastructure, human settlements, and environmental assets.

In this study a mudflow evaluation in a tropical silty soil is presented, from the implementation of a laboratory scale experiment with calibration and validation from numerical models to estimate rheological parameters of the material. Moreover, the calibration data of the numerical model were used for a real case study, simulating the slip flow occurred in Yangbaodi, in the southeast of China, occurred on September 18, 2002. The calibrated model reproduced in an appropriate way the movement of the flow at laboratory scale, and for the aforementioned case study, some differences in the final length of deposition were noticed, achieving interesting results that lead the use of the calibrated model towards the estimation of risks due to the mudflow occurrence.

In Computational Flow Dynamic (CFD) modeling, the equations governing the flow are derived from the theory of multiphase flows. This model originates from the conservation equations of mass, momentum and energy.

2.1. Governing equations

As part of the modeling process, the differential form of the law of conservation of mass and the momentum or linear momentum equations that describe the amount of motion of a body must be used (Versteeg and Malalasekera 2007). Since the turbulent transport process cannot be calculated by an exact method, it must be approximated by a turbulence model that allows the turbulent transport quantities to be related to the mean flow field. In the Reynolds-averaged procedure, the Navier-Stokes equations are averaged over a period of time longer than the largest characteristic period of the flow motion yielding the Reynolds-averaged Navier-Stokes equations (RANS) (Versteeg and Malalasekera 2007; Galván et al. 2011):

$$\frac{\partial {U}_{j}}{\partial {x}_{j}}=0$$

$$\rho \frac{\partial {U}_{i}}{\partial t}+\rho \partial {U}_{j}\frac{\partial {U}_{i}}{\partial {x}_{j}}=-\frac{\partial P}{\partial {x}_{i}}+\rho {g}_{i}+\frac{\partial }{\partial {x}_{j}}({\tau }_{ij}-\rho \stackrel{-}{{u}_{i}{u}_{j}})$$

Where: U_i corresponds to the mean flow velocity vector, x_i, direction in space, t, time, P, flow pressure, g_i, gravitational acceleration, τ_i,j, mean flow deviatoric stress tensor, u_i, low velocity vector for fluctuations and ($-\rho \stackrel{-}{{u}_{i}{u}_{j}})$_, Reynolds stress tensor. The latter cannot be expressed as a function of the mean flow and must be related to known quantities using a turbulence model (Galván et al. 2011; Ahsan 2014).

2.2. Turbulence model

Most flows of high interest in engineering, in this case muds, are turbulent in nature. Therefore, this flow regime is not only of theoretical interest, but also a problematic source for different case studies, where it is necessary to consider the effects of turbulence to better understand the behavior of the flow (Versteeg and Malalasekera 2007; Tu et al. 2013). There are several turbulence models that CFD tools use for their resolution. The best known of all is the standard k-ε model (Zhu and Kuznetsov 2005). Because the complexity of turbulence in most engineering flow problems, if the use of simple formulas is included, it is possible to develop similar transport equations to accommodate the turbulent quantity k and other turbulent quantities, one of which is the turbulent energy dissipation rate ε (Galván et al. 2011; Tu et al. 2013; Trujillo-Vela et al. 2020). Based on the standard model, the RNG (Renormalization Group) and RAE (Realizable) models were developed for k and ε (Shih et al. 1995; Leveque 2004; Galván et al. 2011; Ahsan 2014).

RAE realizable model, which corresponds to the model that will be used in this research, was developed by Shih et al. (1995), is more recent with respect to the standard and RNG models. It positively satisfies the constraints on the related quantities in the calculation of Reynolds stresses, in accordance with the physics of turbulent flows. The model transport equations for k and ε in the realizable modal are:

$$\rho \frac{\partial k}{\partial t}+\rho {U}_{i}\frac{\partial k}{{\partial x}_{i}}=\frac{\partial }{\partial {x}_{j}}\left[\left(\mu +\frac{{\mu }_{t}}{{\sigma }_{k}}\right)\frac{\partial k}{\partial {x}_{j}}\right]+{\tau }_{ij}\frac{{\partial U}_{i}}{\partial {x}_{j}}-\rho \epsilon$$

$$\rho \frac{\partial \epsilon }{\partial t}+\rho {U}_{i}\frac{\partial \epsilon }{{\partial x}_{i}}=\frac{\partial }{\partial {x}_{j}}\left[\left(\mu +\frac{{\mu }_{t}}{{\sigma }_{\epsilon }}\right)\frac{\partial \epsilon }{\partial {x}_{j}}\right]+\rho {C}_{1}{S}_{\epsilon }-\rho {C}_{2}\frac{{\epsilon }^{2}}{k+\sqrt{v\epsilon }}$$

Where: the model constants are C₂ = 1.9, σ_ε = 1.2, σ_k = 1.0, C₁ = max [0.43; (η/η + 5)] (Shih et al. 1995). Therefore, η is given by:

$\eta =\frac{Sk}{\epsilon }$	(5)
$S=\sqrt{2{S}_{ij}{S}_{ij}}$	(6)

Where: S is the modulus of the mean strain rate tensor. The model, unlike the standard one, uses a different formulation for turbulent viscosity, and uses a new model for the dissipation ε (See Eq. 8), where the production term (second on the right) does not involve the production of k. Furthermore, the denominator of the dissipation term (the third one on the right) never vanishes, even if k vanishes or is less than zero. Therefore, these features allow the model to satisfy the mathematical constraints on the normal stresses, according to the physical behavior of turbulent flows (Shih et al. 1995; Galván et al. 2011).

2.3. Multiphase model

To correctly model mud flow, a multiphase air-mud model must be implemented that considers two phenomena: the free-flowing surface or the air-mud interface, which can be realized through a VOF (Volume of fluid) model and a mixing model.

2.3.1. VOF Model

To simulate a multiphase flow, this model uses a unique set of conservation equations for the entire flow domain. This set must consider the differences in the properties of the fluids involved, as well as the surface tension at the interfaces. Eq. 1 does not change, and the motion of the flow is given by the Navier-stokes equations for a compressible fluid presented in Eq. 2 (Jeong and Yang 1998; Lopez 2017). To this, an additional term is added to account for the surface tension present at the air-mud interface (Manninen and Taivassalo 1996):

$$T=\int \sigma {k}^{{\prime }}{n}^{{\prime }{\delta }^{\beta }}\left(x-{x}^{{\prime }}\right)ds{\prime }$$

Where: σ is the surface tension coefficient, k’ is the curvature of the free surface in the 2D domain, n’ the unit vector on the surface δβ is the delta Dirac function, where β = 2 for 2D domains and β = 3 for 3D domains, x is the point on the coordinate axis where the equation is evaluated and x’ is the point on the free surface, such that if x = x’ then the function δ is nonzero, guaranteeing that the term T is evaluated only at the fluid interface, ds’ corresponds to the free surface differential. In equations 1 and 2 ρ is replaced by (Manninen and Taivassalo 1996):

$${\rho }_{m}=\alpha {\rho }_{w}+(1-\alpha ){\rho }_{a}$$

Where: α corresponds to the mud volume fraction, defined as the ratio of the mud volume to the total volume of the mixture, ρ_w is the mud density and ρ_a is the air density. Finally, it is necessary to take into account an additional transport equation for the fluid volume fraction:

$$\frac{\partial \alpha }{\partial t}+\frac{\partial }{\partial {x}_{i}}\left(\alpha {U}_{i}\right)=0$$

2.3.2. Mixture Model

The mixing model uses the conservation of mass and momentum equations for the air-mud mixture, and the continuity equation for the dispersed phase (Manninen and Taivassalo 1996). As the momentum equation is not considered for the dispersed phase, the relative or slip velocity, which is defined as the velocity of the secondary phase p, relative to the velocity of the primary phase q, is considered with the following algebraic expression:

$$\overrightarrow{{v}_{pq}}= \overrightarrow{{v}_{p}}-\overrightarrow{{v}_{q}}$$

The relative velocity can be written as:

$$\overrightarrow{{v}_{pq}}= \frac{{\tau }_{p}}{FD}\frac{({\rho }_{p}-{\rho }_{m})}{{\rho }_{p}}\overrightarrow{\alpha }$$

Where:

$${\tau }_{p}= \frac{{\rho }_{p}{{d}^{2}}_{p}}{{18\mu }_{q}}$$

It is the particle (air bubble of the dispersed phase), FD is the drag function, ρp is the density of the dispersed phase and ρ_m the density of the mixture, α, the particle acceleration of the secondary phase, dp, the particle diameter in the secondary phase and µ_q, the dynamic viscosity of the primary phase. FD is calculated from the Schiller & Newman approach (Manninen and Taivassalo 1996):

$$FD=\left\{\begin{array}{c}1+0.15{Re}^{0.687} Re\le 1000\\ 0.018Re Re>1000\end{array}\right.$$

Where the Reynolds number is given by:

$$Re=\frac{{\rho }_{q}\left|\overrightarrow{{v}_{pq}}\right|{d}_{p}}{{\mu }_{q}}$$

In order to carry out the processes that are fulfilled in the computational flow dynamics (CFD) tools, it is necessary to start with the approximation of a continuous variable in a number of finite points. This corresponds to a spatial discretization process, where, starting from the problem, the equations governing the flow are posed and solved (Leveque 2004; Versteeg and Malalasekera 2007). For this case, the discretization method used was the finite volume method.

2.4. Flow viscosity

According to Leveque (2004); and Versteeg and Malalasekera (2007), the equations governing the flow contain as additional unknowns the viscous stress τ_ij. The most useful forms of the conservation equations for mud flows are obtained by introducing a suitable model for such viscous stresses τ_ij. In many flows, viscous stresses can be expressed as functions of the local strain rate or strain rate.

Since the material behaves as a non-Newtonian fluid, it is essential to consider the rheological behavior of the material as a complement to the development of the model, since the flow of mud, characterized as a highly viscous material, can present variations in its physical properties that influence its progress through any surface (Mahajan and Budhu 2006; Huang & Aode, 2009; Carotenuto et al. 2015; Kameda and Hirauchi 2018). Therefore, the behavior of the material is not linear, but follows a behavior through a mathematical model, which describes the relationship between the shear stress and the strain rate of the material. Studies such as Huang and García (1998); Bartosik (2010); Sawyer et al. (2012); and Han et al. (2019), have used the Herschel-Bulkley model to study the rheological behavior of mudflow to better predict its behavior. This model is given by equation:

$$\tau ={\tau }_{y}+K{\gamma }^{n}$$

Where: τ is the shear stress, τ_y is the yield stress, K is the apparent viscosity, γ is the shear strain rate and n describes the rate of viscosity change. The rheological behavior of the material according to Herschel Bulkley has the form presented in Fig. 1, differing from the behavior of a Newtonian fluid, which presents a linear behavior. Thus, this nonlinear behavior shows that, depending on the type of soil and moisture content, each mud flow presents a different nonlinear behavior.

For the development of this research, laboratory and computational tests were carried out with numerical modeling tools, using the following methodology:

3.1. Soil characterization

For the execution of the tests, a silty soil of residual origin from the northwestern region of Colombia was used, with which the mud was prepared by mixing the soil with different water contents. This material was tested for granulometry (ASTM D422), consistency limits (ASTM D423 and ASTM D424) and specific gravity (ASTM D854). The results are shown in Table 1.

Table 1

Soil characterization
	Value
Liquid limit (%)	43.00
Plastic limit (%)	13.00
Specific gravity	2.71
Particle size: 2 mm (% fine)	-
Sieve N°40: 425 µm (% fine)	-
Sieve N°200: 75 µm (% fine)	84.25
USCS classification (ASTM D2487)	ML

As part of the input data in CFD tools for flow modeling, one of the most important parameters in material characterization is the flow viscosity. Although, muds are considered as a non-Newtonian fluid in this research the mud was also considered as a Newtonian fluid, in order to compare its influence on the flow behavior.

3.1.1. Viscosity as a Newtonian Fluid

Initially, a linear behavior of the flow was considered, where a viscosity value was obtained for each mudflow according to its moisture content. For this purpose, soil muds were prepared with various moisture contents in order to obtain a viscosity curve as a function of moisture, considering that Stokes' law is fulfilled. Mudflows were prepared at 60%, 65%, 70% and 80% moisture. These samples were poured into a test tube. After this, a lead sphere was weighed and released to flow along a length L of the material inside the test tube. (Fig. 2). The elapsed time in the lead ball advance was taken and recorded for the calculation of the viscosity value by means of the following expression:

$$\eta =\frac{2*g*\left({\rho }_{p}-{\rho }_{m}\right)*{r}^{2}}{9*\left(\frac{L}{t}\right)}$$

Where: η corresponds to the viscosity of the material, g corresponds to the value of gravity, ρ_p to the density of the lead sphere, ρ_m to the mud density, r to the radius of the lead sphere, L the length the sphere traverses vertically along the material, and t the elapsed time in which the sphere traverses length L.

3.1.2. Viscosity as a non-Newtonian fluid

To consider the mud as a non-Newtonian material, a regression analysis was performed by a trial-and-error process, obtaining values of the rheological parameters corresponding to the Consistency K, yield stress τ_y and rate of change of viscosity n, corresponding to the model proposed by Herschel-Bulkley. The regression analysis started from data found in the literature for similar soils in the works proposed by Bartosik (2010); Sawyer et al. (2012); Han et al. (2019).

To validate the data, experimental models were made and compared with numerical models. For this purpose, a scaled prototype similar to a Bostwick consistometer was developed, based on the experiments of Balmforth et al. (2007), with dimensions 1 m length, 0.2 m height and 0.2 m width, with slope control at the bottom and a gate at 0.1 m to give a capacity of 0.1 m x 0.2 m x 0.3 m of mud volume before flowing, as shown in Fig. 3. Soil moisture contents of %, 65%, 70% and 80% were used. Finally, the gate was opened to let it slide through the consistometer, using slopes of 5%, 10%, 15% and 20%. The flow runout and mud deposition heights were recorded using the methodology described in Section 3.3.

After, the regression analysis was established for the rheological parameters of the Herschel-Bulkley model corresponding to each soil. Numerical models were made to simulate what happened in the laboratory by varying each parameter. The models made to validate the laboratory observations were configured in a CFD tool under the methodology shown in section 3.3.

The results obtained were used to evaluate the performance of the model according to the methodology proposed by Ritter and Muñoz-Carpena (2013), where a check was performed to quantify the prediction error in terms of the units of the variable calculated by the model. The root mean square error (RMSE) is selected and its definition is given by:

$$RMSE=\sqrt{\frac{{\sum }_{i=1}^{N}{\left({O}_{i}-{P}_{i}\right)}^{2}}{N}}$$

Where O_i and P_i represent the sample (of size N) containing the observations and model estimates, respectively. This value varies between 0 and ∞, where, for an RMSE = 0, it indicates that the model fit is perfect according to the experimental data obtained (Ritter and Muñoz-Carpena 2013). The sample corresponds to the number of deposition heights obtained in the experiment along the entire length. Third, a goodness-of-fit indicator is calculated given by:

$$NSE=1-{\left(\frac{RMSE}{SD}\right)}^{2}$$

Where SD represents the standard deviation of the laboratory observations. Therefore, NSE represents the unit complement of the ratio of the mean square error of the observed values versus the modeled values and the variance of the observations. The efficiency coefficient takes values -∞<NSE < 1 (Ritter and Muñoz-Carpena 2013).

It is then possible to verify the goodness of fit of the model obtained with respect to the observations in the laboratory from the relationship between SD/RMSE and the values obtained for NSE, as shown in Table 2.

Table 2

Criteria for evaluation of model performance, adapted (Ritter and Muñoz-Carpena 2013)
Performance classification	Interpretation	NSE
Very good	[SD/RMSE] > 3.2	> 0.90
Good	2.2 < [SD/RMSE] < 3.2	0.80–0.90
Acceptable	1.2 < [SD/RMSE] < 2.2	0.65–0.80
Unsatisfactory	[SD/RMSE] < 1.7	< 0.65

3.2. Laboratory channel tests

For the laboratory tests, a second scale model was made, consisting of an open channel 3 m length, 0.5 m height and 0.7 m width, with slope control at the bottom and a gate 0.25 m from the start to give a capacity of 0.25 m x 0.5 m x 0.7 m of mud volume to flow after opening the gate (Fig. 4 and Fig. 5).

The tests were carried out with slopes of 5%, 10%, 15% and 20%, considering mud with 60% and 80% moisture. For the tests with 60% moisture mud, the material was placed up to a height of 25 cm, and for the tests with 80% moisture mud, the material was placed up to a height of 40 cm. For the material at 60% moisture, a height of 25 cm was defined in order to guarantee that the channel gate would open properly, since higher heights with this type of material generated difficulties with the opening of the gate, which is why it was defined.

To carry out the measurement with a laser distance meter, a grid system was established in the upper part of the channel, using 10 cm x 10 cm cells (Fig. 6). through of this system, a distance meter was located in the center of the cell to verify the deposition heights in that area, then, considering the number of cells in the flow runout, it was possible to determine the final length on the channel surface. The deposition height was calculated as:

$$HL=H-n$$

3.3. Computational dynamic modeling of mudflow

The mudflow from the channel was calculated using a student version CFD package, in this case FLUENT 2020 R2. A two-dimensional virtual geometric model of the channel was implemented and a numerical model was assigned to simulate the flow behavior. Initially, a calibration of the numerical model was performed, for which the numerical results obtained in the software were compared with the results measured in the laboratory scale model for the 5% slope. Subsequently, the model was validated by comparing the test results with 10%, 15% and 20% slope. For the construction of the 2D model, the flow profile passing through section D was simulated according to the grid system generated in Fig. 6 during laboratory tests for measurement.

3.3.1. Geometric model and boundary conditions

The geometry of the model is shown in Fig. 7 with the initial boundary conditions established. Thus, three boundaries and two zones were established: Wall boundary where the mud rests and flows after the gate opening, atmosphere boundary in the upper part, and outlet boundary in the right part of the geometry, which is where the mudflow exits the channel. Additionally, an initial condition is established, with a zone corresponding to the mud at rest before flowing and another zone that is occupied by air.

For the atmosphere and outlet type boundaries, a pressure outlet boundary condition is used. A gauge atmospheric pressure of 0 Pa is set and a backflow air volume fraction of 1.0 is defined to ensure that in cases of backflow, only air enters the domain (Lopez 2017). For the wall-type boundary, a stationary wall boundary condition is set at the bottom of the channel with a no-slip condition. A wall roughness height ks = 0 mm is specified, considering that for the physical model one has a smooth bottom in which the roughness is negligible.

As shown in Fig. 8, 2D geometry of the model is composed of a zone where the material with volume fraction = 1 and velocity = 0 is at rest. The moment at which the mud starts to flow along the surface corresponds to the moment at which the channel gate opens at laboratory scale.

3.3.2. Mesh generation for the model

To run the numerical model in the CFD tool according to the 2D geometry of the channel, meshes were made for each case of slope and material properties, to be solved by the finite volume methodology. The characteristics of the developed meshes are presented in Table 3. These features are the final results of a mesh sensitivity analysis, gradually increasing their resolution until the model results were shown to be independent of cell size. This procedure was done for the 5% mesh slope as part of the calibration process. Finally, for the area where the material is found before flowing along the surface, a refinement of the cells to 1x10^− 3 m was made with respect to the rest of the domain, which is where the air is found, where cell sizes of 2 x10^− 3 m are present, as shown in Fig. 9. Considering this, it can be guaranteed that the material properties will not be affected by the resolution of the mesh in the air domain, generating a gradual transition in the cell size, which allows for better flow conditions when running the model.

Table 3

Number of elements per mesh according to slope and moisture content of the material
Slope	Soil moisture	# Elements
5%	60%	62970
10%	60%	68310
15%	60%	73140
20%	60%	78330
5%	80%	65830
10%	80%	71130
15%	80%	76810
20%	80%	80990

A refinement was also made over the entire length of the bottom of the channel, which is the area where the material will flow, also with the objective that this area, which initially has air and will be occupied by the material, does not influence its properties. For this purpose, an inflation algorithm was applied with 15 layers and a growth rate of 1.2, with a height for the first layer of 5x10^− 4 m, in order to cover the boundary layer (ANSYS 2013).

3.3.3. Mesh generation for the model

The mixture model was used as the default model for the development of the models. Figure 9 shows how the volume fraction occupies the air present in the initial condition and a transition of the flow into the air is observed, thus, in order for everything to happen correctly, the properties of the model shown in the Table 4.

Table 4

Model characteristics
Initial setting	Flow type	Transient
Initial setting	Gravity (m/s²)	(0, -9.81)
Multi-phase flow	Primary phase (Mud)	Density: (Table 8)
Multi-phase flow	Secondary phase (Air)	Particle size: 5 x 10^− 3 m
Turbulence model	Model k – ε realizable (RAE)
Viscosity model	Herschel – Bulkley parameters (Table 8)
Solution method	Coupling velocity-pressure	SIMPLE
	Spatial discretization (Pressure)	PRESTO!
	Momentum	QUICK
	Volumetric Fraction
	Kinetic turbulent energy
	Specific dissipation rate
Stopping criteria	Flow velocity = 0.0 m/s

The results obtained in the modeling process will be initially compared with the results obtained in the laboratory to observe how well the model is reproducing the flow runout and the deposition height along the channel surface. Also, the models under the non-Newtonian flow and Newtonian flow conditions will be compared to observe what differences there are between them and what the scope of each of them can be when required to model real case studies involving mudflows.

3.4. Case Study

To validate all the results obtained from the models corresponding to the channel, a numerical modeling process was carried out on a full-scale material flow. The case study corresponding to the catastrophic rainfall-induced flowslide in Yangbaodi (Shenzhen), in the southeast of China, occurred at 1:40 p.m. on September 18, 2002, after the slope experienced prolonged rainfalls. About 25,000 m³ of fill soils consisted of completely and intensely decomposed granite, traveling for a distance of about 140 m at a considerably high speed, which was Li et al., (2012). using the Discrete Element Particle Method (DEM).

The initial geometry for the 2D modeling that was used in the works by Li et al. (2012) and Zhang et al. (2015) is presented in Fig. 10. It is important to note that the sliding material is assumed to behave as a continuous medium, as considered in the laboratory scale tests. Thus, the model configuration parameters in ANSYS FLUENT are the same as those described in Section 3.4, and the material properties are described in Table 5.

Parameters were used to simulate the behavior of the soil under the model proposed by Herschel-Bulkley through data found in the literature that simulate the rheological behavior of the material. These parameters were determined according to the work carried out by (Jeong 2019), which had very similar characteristics to the materials present in the Yangbaodi landslide.

Table 5

Material properties for full-scale model
Parameter	Value
Density (kg/m³) (Zhang et al. 2015)	1,400
Exponent Herschel-Bulkley n (Jeong 2019)	1.14
Consistency K (Jeong 2019)	0.83
Fluency Stress τ_y (Pa) (Jeong 2019)	23.00

After performing the rheological characterization tests of the material, the laboratory tests in the channel and the modeling processes described above, the following results and comparisons were obtained:

4.1. Viscous behavior of the material considering the material as a Newtonian fluid.

Following the procedure described in section 3.1.1, tests were carried out to measure the viscosity of the mixtures, considering that the analysis of the materials also considered the mudflow as a Newtonian. Figure 11 presents the results of viscosity and density for different moisture contents.

The results show a decrease in the viscosity value when the water content in the material increases because a higher amount of water decreases the attraction between particles, which decreases the shear force present. This decreases the resistance of the material to flow and, therefore, implies greater runout along a surface. When the material reaches humidities higher than 70%, the viscosity behavior tends to be asymptotic.

With the density of the material the behavior is similar, since there is a decrease as the water content increases, which contributes less weight in a control volume. Therefore, the density is directly proportional to the viscosity value, where higher densities imply a higher resistance to flow and therefore higher viscosity values.

4.2. Viscous behavior of the material considering the material as a non-Newtonian fluid.

According to the procedure described in section 3.1.2, viscosity values of the material were determined considering it as a non-Newtonian fluid. The values corresponding to the consistency K, yield stress τ_y and rate of change of viscosity n, of the Herschel-Bulkley model were calibrated through the regression analysis and are shown in Table 6 according to the material moisture content.

Table 6

Rheological parameters for the material according to the Herschel-Bulkley model
Moisture content	Fluency stress τ_y (Pa)	Consistency K	Exponent Herschel-Bulkley n
60%	100	80	0.45
65%	72	50	0.45
70%	21.5	10	0.46
80%	15.5	5	0.46

According to Table 3 and Fig. 16, When the moisture content increases, the interaction between adjacent particles decreases and, therefore, the yield stress decreases. Therefore, knowledge of composition or water content alone is insufficient to predict yield strength and mud flow behavior. As material rich in silt content is shown, there could be changes in the flow velocity that could be higher or lower, which directly depends on the grain size, so constant Herschel-Bulkley exponent (n) values can be noted, since in this research the changes are linked to the conditions of the experiment under a single material. This agrees with that presented by (Sawyer et al. 2012).

Models results with the rheological data of the calibrated material are shown in Fig. 15, where the deposition heights and flow reach length for each consistometer slope condition and material moisture content are described. Additionally, it is compared with the observations obtained in the laboratory.

The iterative process was carried out until the models were adequately adjusted to the results of the consistometer tests and a good efficiency was obtained, in accordance with the criteria of Table 2. With equations 17 and 18, the SD/RMSE ratio and the NSE value were calculated, and compared with the parameters of Table 3 for each of the models proposed, and the results are shown en in Table 7.

Figure 12 indicates a good reproduction of the model over the distance traveled and deposition heights obtained, the model performance results show that the model satisfactorily matches the observations obtained in the test, indicating that the Herschel-Bulkley model parameters used satisfactorily reproduce the flow movement along the surface.

Table 7

Performance results of the model used for the consistometer
Test	Slope	Soil moisture	SD/RMSE	NSE	Classification
1	5%	60%	11.58	0.99	Very good
2	5%	65%	8.23	0.99	Very good
3	5%	70%	3.73	0.93	Very good
4	5%	80%	2.86	0.88	Good
5	10%	60%	4.74	0.96	Very good
6	10%	65%	8.97	0.99	Very good
7	10%	70%	4.14	0.86	Very good
8	10%	80%	4.65	0.90	Very good
9	15%	60%	5.39	0.97	Very good
10	15%	65%	7.42	0.98	Very good
11	15%	70%	4.14	0.94	Very good
12	15%	80%	4.65	0.88	Very good
13	20%	60%	2.50	0.84	Good
14	20%	65%	2.24	0.80	Good
15	20%	70%	2.82	0.87	Good
16	20%	80%	4.65	0.95	Very good

Figure 13 shows the results obtained from the regression analysis performed, which describe the rheological behavior of the materials, considering moisture contents of 60%, 65%, 70% and 80%. On the other hand, the model was able to determine the change in shear stress as a function of shear rate for each moisture content, whose variation is represented on the abscissa axis. For each material, according to its water content, the steady-state yield stress can be determined for each strain rate, defining τy at zero strain rate. The results were fitted with the model proposed by Herschel-Bulkley.

The material with 60% or 65% moisture content presents a greater increase in shear strength as the deformation rate increases, while in materials with 70% and 80%, a behavior with a lower growth rate is observed due to the amount of water they have, which decreases the internal shear forces and allows a greater fluidity of the material, where, from 70%, the material flows more easily along a surface and allows the solids content to be dragged easily by the volume of water present.

4.3. Channel tests and numerical models

Tests were carried out in the channel as described in section 3, where it was found that the material for 60% moisture content, regardless of the slope used in the channel, stops between 0.88 and 1.2 m in length along the surface. For 80% moisture, the flow runout to the end of the channel and an accumulation of material are generated at the end of the channel, especially in the last 0.2 m of the channel.

Models Results, according to slope, material moisture, fluid density and flow type, are described in the following Table 8.

Table 8

Material characteristics for the models to be made according to slope and moisture content
Model	Slope	Soil moisture	Flow type	Density(kg/m3)	Viscosity Pa*s	Fluency stress (Pa)	Consistency K	Exponent Herschel Bulkley n
1	5%	60%	Non Newtonian	1487.20	N/A	100.0	80	0.45
2	5%	80%	Non Newtonian	1424.45	N/A	15.5	5	0.46
3	5%	60%	Newtonian	1487.20	90.23	N/A	N/A	N/A
4	5%	80%	Newtonian	1424.45	0.887	N/A	N/A	N/A
5	10%	60%	Non Newtonian	1487.20	N/A	100.0	80	0.45
6	10%	80%	Non Newtonian	1424.45	N/A	15.5	5	0.46
7	10%	60%	Newtonian	1487.20	90.23	N/A	N/A	N/A
8	10%	80%	Newtonian	1424.45	0.887	N/A	N/A	N/A
9	15%	60%	Non Newtonian	1487.20	N/A	100.0	80	0.45
10	15%	80%	Non Newtonian	1424.45	N/A	15.5	5	0.46
11	15%	60%	Newtonian	1487.20	90.23	N/A	N/A	N/A
12	15%	80%	Newtonian	1424.45	0.887	N/A	N/A	N/A
13	20%	60%	Non Newtonian	1487.20	N/A	100.0	80	0.45
14	20%	80%	Non Newtonian	1424.45	N/A	15.5	5	0.46
15	20%	60%	Newtonian	1487.20	90.23	N/A	N/A	N/A
16	20%	80%	Newtonian	1424.45	0.887	N/A	N/A	N/A

4.3.1. Numerical model results considering the material as a Newtonian fluid

Different scenarios were numerically simulated and plotted. Figure 14 shows the runout and the deposition height of the flow obtained in the numerical simulation considering all the slopes described in Table 8 and that the flow is Newtonian with moisture contents of 60% and 80%. In addition, the deposition heights and flow runout obtained in the laboratory tests are presented as dotted lines.

The results obtained in the numerical simulations tend to reproduce what occurs in the physical model. However, when comparing each test with their respective simulations, some variations are found during the flow path, especially when analyzed for 60% moisture, where the deposition length exceeds between 20 cm and 30 cm the distance that the flow traveled in the laboratory. This difference tends to be the same for each slope used, showing that for higher slopes and material with 60% moisture there is a tendency to flow a greater distance, which occurs due to the action of gravity. For 80% moisture, the numerical simulation yields a flow profile with the same range as in the laboratory, with an average difference in the height of deposition along the surface of approximately 2 cm above the flow profile obtained in the laboratory tests.

4.3.2. Numerical model results considering the material as a non-Newtonian fluid

Figure 15 shows the results obtained in the simulations carried out to reproduce again the phenomenon obtained in the laboratory, but this time modeling the material as a non-Newtonian fluid, considering the calibrated viscosity parameters described in the Table 8 and all slopes. The flow profiles obtained in these simulations show that the final deposition heights, as well as the final runout agree with the results obtained at laboratory level. When analyzing the flow profile in all slopes, when the material has a water content of 60%, it can be observed that the difference in the final length of deposition varies between 5cm and 10cm more than the displacement of the material in the physical model.

Figure 15c and Fig. 15d show that this difference is greater with respect to the results shown in Fig. 15a and Fig. 15b, which also coincides with the increase in the slope of the channel. However, it could be said that the slope, when the material presented 60% moisture, and according to the deposition lengths obtained, was not influential on the flow advance, which, after stopping, remained between 0.9 m and 1.2 m for the laboratory tests and between 0.9 m and 1.3 m for the results of the numerical simulations. Thus, it would be important to verify larger and variable slopes, based on real case studies, as shown in Section 4.4, in order to increase the validity of the results obtained.

For materials with 80% moisture, the results show graphically that there is a tendency to reproduce not only the flow motion along the surface, but also the final length and final deposition heights, which extend along the channel. This means that the rheological parameters of the material have an influence on the flow motion and largely control how they move along the surface.

4.3.3. Model Analysis

After qualitatively analyzing the results, Equations 17 and 18 were applied, the SD/RMSE ratios and NSE values are presented to verify the goodness of fit with the laboratory tests.

Table 9

SD/RMSE ratio for goodness of fit of the model
Model	Slope	Soil moisture	Flow type	SD/RMSE	NSE	Classification
1	5%	60%	Non Newtonian	2.51	0.84	Good
2	5%	80%	Non Newtonian	1.27	0.38	Acceptable
3	5%	60%	Newtonian	1.74	0.67	Acceptable
4	5%	80%	Newtonian	0.63	-1.54	Unsatisfactory
5	10%	60%	Non Newtonian	2.96	0.89	Good
6	10%	80%	Non Newtonian	1.34	0.44	Acceptable
7	10%	60%	Newtonian	1.85	0.71	Acceptable
8	10%	80%	Newtonian	0.63	0.39	Unsatisfactory
9	15%	60%	Non Newtonian	3.30	0.91	Very good
10	15%	80%	Non Newtonian	1.66	0.64	Acceptable
11	15%	60%	Newtonian	2.64	0.86	Good
12	15%	80%	Newtonian	1.62	0.62	Unsatisfactory
13	20%	60%	Non Newtonian	1.73	0.67	Acceptable
14	20%	80%	Non Newtonian	1.97	0.74	Acceptable
15	20%	60%	Newtonian	1.08	0.14	Unsatisfactory
16	20%	80%	Newtonian	1.64	0.63	Unsatisfactory

According to the results shown, the Table 9, it is observed that when the flow is modeled considering the material as a Newtonian fluid, the efficiency and goodness of fit of the model to the experimental observations tend to be unsatisfactory in most models, especially for 80% moisture. On the other hand, models based on non-Newtonian fluids tend to reproduce what is observed in the laboratory, confirming the theory that mudflow is a non-Newtonian fluid. However, it is important to consider that due to the conditions of the test in the channel for mud with 80% moisture, there was a loss of material at the end of its course, due to the small variation in the shear stress due to its water content, which caused faster movements and stronger shocks at the end of the channel, limiting some measurements in the deposition heights.

The model for non-Newtonian flows reproduces the models better than the model for Newtonian flows, probably because by simulating with calibrated rheological parameters, the CFD tool is able to assimilate a viscous material behavior according to the variation in its consistency and shear force, improving the detail in the reproduction of its motion.

4.4. Model results in case study: Yangbaodi Landslide

The case study was simulated using the configuration of the models that were used to reproduce the laboratory tests, considering the behavior of the material as that of a non-Newtonian viscous fluid. The boundary conditions and model configuration were the same as those used in section 3.3.

Depending on the material and model configurations, a flow was produced from the initial position of the material. This flowed along the entire surface, achieving the breakthrough length shown in the Fig. 16, presenting at the abscissa 0 + 175 of the flow a height of 135 m with respect to the abscissa 0 + 327, where the flow presents a height of 122 m, which is where all the material was deposited longitudinally when this phenomenon occurred. Analyzing the final distance and deposition heights, the model results show a greater accumulation of material in the first 90 m, after which the deposition heights are reduced until reaching the final point where the flow stopped at the front.

The total deposition length of the flow according to the numerical simulation was 146.56 m, while the deposition length for the flow at Yangbaodi was 123.43 m. The model shows differences at the 0 + 190 abscissa, stopping at approximately the 0 + 179.5 abscissa, and distributing the material to a distance 9.9 m more than the actual deposition distance of the flow. Although these differences are present, the model adjusts and reproduces the movement of the material generating a distance traveled. The variations probably occur because properties of a material with similar characteristics were used, but which may have differences in their physical conditions, considering that the parameters K and τ_y of the Herschel-Bulkley model define the rate of variation in the shear strength, which could have generated these variations in the flow path and therefore the way it was deposited. (Fig. 17).

According to the results obtained, it is evident that CFD is a useful tool for the modeling of mudflows or flow-type landslides, which could be used in the evaluation of hazard and risk due to landslides. By means of these techniques, travel distances, deposition height, and variables such as pressures on structures, drag forces, among others, can be determined.

The "k-ε realizable" turbulence model adequately simulates the mudflow with respect to the results obtained in the experimental tests. When considering the material as a non-Newtonian fluid or as a Newtonian fluid within the model, some variations are presented in some coordinates of the flow contrasted with the final deposition heights and runout where it is observed that according to the SD/RMSE ratio and the NSE value, for models with non-Newtonian flows the experimental results are captured in a better way. Thus, for case studies, where there are important variations in the surface where the flow advances, the difference in the results given by both models could increase, so that the turbulence model considering the mudflow as a Newtonian fluid could fail to adequately predict these processes.

The slope of the terrain plays a crucial role in determining the extent of mudflow, including the height at which the mud ultimately settles. As the slope becomes steeper, the mudflow can travel greater distances as it is influenced by the force of gravity. Consequently, when the mudflow covers longer distances, the final deposition heights tend to be lower because the material is spread out over a larger area along its path. This effect is especially pronounced in mudflows with a high moisture content of around 80%, as the increased water content reduces the material's viscosity and the shear stress between individual particles. This reduction in consistency allows the mudflow to advance over longer distances before coming to rest.

The used model adequately reproduces the flow movement over an irregular surface. In this case, when simulating the event ocurred in the Yangbaodi, the flow presented a final deposition on the area where the phenomenon occurred. The variations obtained occurred due to the variation in the data governing its rheological behavior, so it is necessary for the simulation of case studies to extend the knowledge of the material, in order to have a better perception of the final results of the models.

In this research, flow velocity and pressure were also modeled. With the velocity variable, it was possible to validate that the flow effectively stopped at a certain point on the surface according to the modeling time. Therefore, it is important for future research to integrate the analysis with variables such as flow pressure to quantify the impact on existing structures within the flow path according to their velocity. This allows to improve the estimation of risks due to the mudflows occurrence.

It is recommended to perform experiments and modeling with mud compositions of other types of soils, which may include more granular materials and materials of lower plasticity, and to increase the range of moisture evaluation, especially in scale tests in the laboratory prototype, so that it can be validated with other materials. Additionally, it is recommended to use longer lengths at the experimental and numerical level, in order to have a greater coverage in the evaluation of the flow forward length. Since with longer path lengths the energy of the sliding flow is dissipated, and it is possible to perceive how the final deposition heights are distributed along its advancing length.

Summarizing, based on the provided results, here are five major outcomes that can be drawn:

Regarding the viscosity and density of the material:

The viscosity of the material decreases as the water content increases, indicating a decrease in resistance to flow.
Higher densities correspond to higher viscosity values, implying greater resistance to flow.
Material with humidities higher than 70% exhibits asymptotic viscosity behavior.

Regarding the viscous behavior of the material as a non-Newtonian fluid:

The yield stress of the material decreases with an increase in moisture content.
Changes in flow velocity depend on the grain size and can vary for materials with different silt content.
The Herschel-Bulkley model parameters used satisfactorily reproduce the flow movement along the surface.

Regarding the channel tests and numerical models:

For material with 60% moisture content, regardless of the slope used, the flow stops between 0.88 and 1.20 m in length along the surface.
Material with 80% moisture content advances to the end of the channel, accumulating material at the end.

Regarding the numerical model results considering the material as a Newtonian fluid:

Simulations tend to reproduce the behavior observed in the physical model, but with some variations in the flow path.
The distance of deposition for material with 60% moisture content exceeds the laboratory measurements, while for 80% moisture content, the simulation matches the laboratory results reasonably well.

Regarding the numerical model results considering the material as a non-Newtonian fluid:

Simulations modeling the material as a non-Newtonian fluid accurately reproduce the flow motion, deposition heights, and runout observed in the laboratory tests.

Acknowledgments. Research data was supplied by research program “Vulnerability, resilience and risk of communities and supplying basins affected by landslides and avalanches”, code 1118-852-71251, project “Functions for vulnerability assessment due to water shortages by landslides and avalanches: micro-basins of southwest Antioquia”, contract 80740-492-2020 held between Fiduprevisora and the Universidad de Medellín, with resources from the National Financing Fund for science, technology, and innovation, “Francisco José de Caldas”.

Conflict of interest. The authors declare that they have no confict of interest.

Abraham MT, Satyam N, Pradhan B, Tian H (2022) Debris flow simulation 2D (DFS 2D): Numerical modelling of debris flows and calibration of friction parameters. Journal of Rock Mechanics and Geotechnical Engineering 14:1747–1760. https://doi.org/10.1016/j.jrmge.2022.01.004
Ahsan M (2014) Numerical analysis of friction factor for a fully developed turbulent flow using k–ε turbulence model with enhanced wall treatment. Beni Suef Univ J Basic Appl Sci 3:269–277. https://doi.org/10.1016/j.bjbas.2014.12.001
ANSYS I (2013) ANSYS Fluent Theory Guide (Release 15.0).
Baggio T, Mergili M, D’Agostino V (2021) Advances in the simulation of debris flow erosion: The case study of the Rio Gere (Italy) event of the 4th August 2017. Geomorphology 381:. https://doi.org/10.1016/j.geomorph.2021.107664
Balmforth NJ, Craster R V., Perona P, et al (2007) Viscoplastic dam breaks and the Bostwick consistometer. J Nonnewton Fluid Mech 142:63–78. https://doi.org/10.1016/j.jnnfm.2006.06.005
Bartosik A (2010) Application of rheological models in prediction of turbulent slurry flow. Flow Turbul Combust 84:277–293. https://doi.org/10.1007/s10494-009-9234-y
Buiskikh AA (2015) Placer-forming flows and mudflows as viscous fluids in river valleys. Russian Geology and Geophysics 56:1749–1756. https://doi.org/10.1016/j.rgg.2015.11.007
Cabral V, Reis F, Veloso V, et al (2023) A multi-step hazard assessment for debris-flow prone areas influenced by hydroclimatic events. Eng Geol 313:. https://doi.org/10.1016/j.enggeo.2022.106961
Carotenuto C, Merola MC, Álvarez-Romero M, et al (2015) Rheology of natural slurries involved in a rapid mudflow with different soil organic carbon content. Colloids Surf A Physicochem Eng Asp 466:57–65. https://doi.org/10.1016/j.colsurfa.2014.10.037
Cascini L, Cuomo S, Pastor M, Sorbino G (2010) Modeling of Rainfall-Induced Shallow Landslides of the Flow-Type. Journal of Geotechnical and Geoenvironmental Engineering 136:85–98. https://doi.org/10.1061/(asce)gt.1943-5606.0000182
Feng SJ, Gao HY, Gao L, et al (2019) Numerical modeling of interactions between a flow slide and buildings considering the destruction process. Landslides 16:1903–1919. https://doi.org/10.1007/s10346-019-01220-9
Galván S, Reggio M, Guibault F (2011) Assessment study of k-ε turbulence models and near-wall modeling for steady state swirling flow analysis in draft tube using fluent. Engineering Applications of Computational Fluid Mechanics 5:459–478. https://doi.org/10.1080/19942060.2011.11015386
George DL, Iverson RM (2014) A depth-averaged debris-flow model that includes the effects of evolving dilatancy. II. Numerical predictions and experimental tests. Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences 470:. https://doi.org/10.1098/rspa.2013.0820
Han Z, Su B, Li Y, et al (2019) Numerical simulation of debris-flow behavior based on the SPH method incorporating the Herschel-Bulkley-Papanastasiou rheology model. Eng Geol 255:26–36. https://doi.org/10.1016/j.enggeo.2019.04.013
Hidalgo C, Vega J (2021) Probabilistic landslide risk assessment in water supply basins: La Liboriana River Basin (Salgar-Colombia). Natural Hazards 109:273–301
Huang X, García M (1998) A Herschel-Bulkley model for mud flow down a slope. J Fluid Mech 374:305–333. https://doi.org/10.1017/S0022112098002845
Huang Z, Aode H (2009) A laboratory study of rheological properties of mudflows in Hangzhou Bay, China. International Journal of Sediment Research 24:410–424. https://doi.org/10.1016/S1001-6279(10)60014-5
Hungr O (1995) A model for the runout analysis of rapid flow slides, debris flows, and avalanches. Canadian Geotechnical Journal 32:610–623. https://doi.org/10.1139/t95-063
Hungr O, Evans SG, Bovis MJ, Hutchinson JN (2001) A review of the classification of landslides of the flow type. Environmental and Engineering Geoscience 7:221–238
Iverson RM (1997) The physics of debris flow, Geophys Journal 35:245–296. https://doi.org/10.1029/97RG00426
Jeong SW (2019) Shear rate-dependent rheological properties of mine tailings: Determination of dynamic and static yield stresses. Applied Sciences (Switzerland) 9:. https://doi.org/10.3390/app9224744
Jeong, Yang (1998) Finite element analysis of transient fluid flow with free surface using VOF (volume-of-fluid) method and adaptive grid. Int J Numer Methods Fluids 26:1127–1154. https://doi.org/10.1002/(SICI)1097-0363(19980615)26:10<1127::AID-FLD644>3.0.CO;2-Q
Kameda J, Hirauchi K ichi (2018) Rheological properties of composite serpentine-brucite suspensions: Implications for mudflow behavior on forearc seamounts. Mar Geol 403:191–196. https://doi.org/10.1016/j.margeo.2018.06.004
Lacerda W, Ehrlich M, Fontoura S, et al (2004) Post-failure mechanics of landslides – flowslides in pyroclastic soils. Landslides: Evaluation and Stabilization/Glissement de Terrain: Evaluation et Stabilisation, Set of 2 Volumes 133:1343–1354. https://doi.org/10.1201/b16816-192
Lee SHH, Widjaja B (2013) Phase concept for mudflow based on the influence of viscosity. Soils and Foundations 53:77–90. https://doi.org/10.1016/j.sandf.2012.12.005
Leveque (2004) Finite Volume Methods for Hyperbolic Problems. Meccanica 39. https://doi.org/10.1023/A:1026256231021
Li P, Shen W, Hou X, Li T (2019) Numerical simulation of the propagation process of a rapid flow-like landslide considering bed entrainment: A case study. Eng Geol 263:105287. https://doi.org/10.1016/j.enggeo.2019.105287
Li WC, Li HJ, Dai FC, Lee LM (2012) Discrete element modeling of a rainfall-induced flowslide. Eng Geol 149–150:22–34. https://doi.org/10.1016/j.enggeo.2012.08.006
Li X, Zhao J (2018) A unified CFD-DEM approach for modeling of debris flow impacts on flexible barriers. Int J Numer Anal Methods Geomech 42:1643–1670. https://doi.org/10.1002/nag.2806
Liu C, Sun Q, Zhou GGD (2018a) Coupling of material point method and discrete element method for granular flows impacting simulations. Int J Numer Methods Eng 115:172–188. https://doi.org/10.1002/nme.5800
Liu W, Yan S, He S (2018b) Landslide damage incurred to buildings: A case study of Shenzhen landslide. Eng Geol 247:69–83. https://doi.org/10.1016/j.enggeo.2018.10.025
López de Mesa Aguilar D (2017) Study of skimming flow in flat slope stepped channels. [online] Universidad Nacional de Colombia Sede Medellín Facultad de Minas Escuela de Geociencias y Medio Ambiente.
Mahajan SP, Budhu M (2006) Viscous effects on penetrating shafts in clays. Acta Geotech 1:157–165. https://doi.org/10.1007/s11440-006-0014-8
Manninen M, Taivassalo V, Kallio S (1996) On the Mixture Model for Multiphase Flow, Technical Research Center of Finland, Espoo, VIT Publications 288.
McDougall S (2017) 2014 canadian geotechnical colloquium: Landslide runout analysis — current practice and challenges. Canadian Geotechnical Journal 54:605–620. https://doi.org/10.1139/cgj-2016-0104
Ortiz-Giraldo L, Botero BA, Vega J (2023) An integral assessment of landslide dams generated by the occurrence of rainfall-induced landslide and debris flow hazard chain. Front. Earth Sci. 11, 1–19. https://doi.org/10.3389/feart.2023.1157881
Pastor M, Tayyebi SM, Hernandez A, et al (2023) A new two-layer two-phase depth-integrated SPH model implementing dewatering: Application to debris flows. Comput Geotech 153:. https://doi.org/10.1016/j.compgeo.2022.105099
Peruzzetto M, Levy C, Thiery Y, et al (2022) Simplified simulation of rock avalanches and subsequent debris flows with a single thin-layer model: Application to the Prêcheur river (Martinique, Lesser Antilles). Eng Geol 296:. https://doi.org/10.1016/j.enggeo.2021.106457
Peruzzetto M, Mangeney A, Grandjean G, et al (2020) Operational estimation of landslide runout: Comparison of empirical and numerical methods. Geosciences (Switzerland) 10:1–35. https://doi.org/10.3390/geosciences10110424
Qingyun Z, Mingxin Z, Dan H (2022) Numerical simulation of impact and entrainment behaviors of debris flow by using SPH–DEM–FEM coupling method. Open Geosciences 14:1020–1047. https://doi.org/10.1515/geo-2022-0407
Ritter A, Muñoz-Carpena R (2013) Performance evaluation of hydrological models: Statistical significance for reducing subjectivity in goodness-of-fit assessments. J Hydrol (Amst) 480:33–45. https://doi.org/10.1016/j.jhydrol.2012.12.004
Sawyer DE, Flemings PB, Buttles J, Mohrig D (2012) Mudflow transport behavior and deposit morphology: Role of shear stress to yield strength ratio in subaqueous experiments. Mar Geol 307–310:28–39. https://doi.org/10.1016/j.margeo.2012.01.009
Shih TH, Liou W, Shabbir A, et al (1995) Transitional dynamics of freely falling discs. Notes on Numerical Fluid Mechanics and Multidisciplinary Design 133:105–116. https://doi.org/10.1007/978-3-319-27386-0_7
Trewhela T, Ihle C, Tamburrino A (2014) Numerical simulations of comminution slurries over complex topographies: Putting together CFD and pipeline integrity. Miner Eng 63:139–148. https://doi.org/10.1016/j.mineng.2014.03.005
Trujillo-Vela MG, Galindo-Torres SA, Zhang X, et al (2020) Smooth particle hydrodynamics and discrete element method coupling scheme for the simulation of debris flows. Comput Geotech 125:103669. https://doi.org/10.1016/j.compgeo.2020.103669
Tu J, Yeoh G, Liu C, Tao Y (2013) Computational fluid dynamics: a practical approach.
Vargas-Cuervo G, Rotigliano E, Conoscenti C (2019) Prediction of debris-avalanches and -flows triggered by a tropical storm by using a stochastic approach: An application to the events occurred in Mocoa (Colombia) on 1 April 2017. Geomorphology 339:31–43. https://doi.org/10.1016/j.geomorph.2019.04.023
Vega JA, Hidalgo CA (2016) Quantitative risk assessment of landslides triggered by earthquakes and rainfall based on direct costs of urban buildings. Geomorphology 273:217–235. https://doi.org/10.1016/j.geomorph.2016.07.032
Vega JA, Hidalgo CA (2017) Advancing Culture of Living with Landslides. Advancing Culture of Living with Landslides. https://doi.org/10.1007/978-3-319-53498-5
Vega, J., Hidalgo, C. (2023) Comparison study of a landslide ‑ event hazard mapping using a multi ‑ approach of fuzzy logic, TRIGRS model, and support vector machine in a data ‑ scarce Andes Mountain region. Arab. J. Geosci. 16:527, 1–26. https://doi.org/10.1007/s12517-023-11627-3
Versteeg, HK, Malalasekera W (2007). An Introduction to Computational Fluid Dynamics: The Finite Volume Method. Pearson Education, Harlow
Wang T, Huang T, Shen P, et al (2023) The mechanisms of high mobility of a glacial debris flow using the Pudasaini-Mergili multi-phase modeling. Eng Geol 322:. https://doi.org/10.1016/j.enggeo.2023.107186
Widjaja B, Lee SHH (2013) Flow box test for viscosity of soil in plastic and viscous liquid states. Soils and Foundations 53:35–46. https://doi.org/10.1016/j.sandf.2012.12.002
Widjaja B, Parahyangan UK (2015) Reccomendation of Viscosity Values for Mudflow. https://doi.org/10.13140/RG.2.1.2690.4169
Xia X, Liang Q (2018) A new depth-averaged model for flow-like landslides over complex terrains with curvatures and steep slopes. Eng Geol 234:174–191. https://doi.org/10.1016/j.enggeo.2018.01.011
Xu J, Huhe A (2016) Rheological study of mudflows at Lianyungang in China. International Journal of Sediment Research 31:71–78. https://doi.org/10.1016/j.ijsrc.2014.06.002
Zhang X, Krabbenhoft K, Sheng D, Li W (2015) Numerical simulation of a flow-like landslide using the particle finite element method. Comput Mech 55:167–177. https://doi.org/10.1007/s00466-014-1088-z
Zheng HC, Shi ZM, Peng M, Yu SB (2018) Coupled CFD-DEM model for the direct numerical simulation of sediment bed erosion by viscous shear flow. Eng Geol 245:309–321. https://doi.org/10.1016/j.enggeo.2018.09.003
Zhu J, Kuznetsov A V. (2005) Forced convection in a composite parallel plate channel: Modeling the effect of interface roughness and turbulence utilizing a k-E model. International Communications in Heat and Mass Transfer 32:10–18. https://doi.org/10.1016/j.icheatmasstransfer.2004.05.019
Zhuang J qi, Peng J bing (2014) A coupled slope cutting—a prolonged rainfall-induced loess landslide: a 17 October 2011 case study. Bulletin of Engineering Geology and the Environment 73:997–1011. https://doi.org/10.1007/s10064-014-0645-1

Download PDF

Editorial decision: Major revisions
30 Sep, 2024
Reviewers agreed at journal
25 Sep, 2023
Reviewers invited by journal
25 Sep, 2023
Editor assigned by journal
13 Sep, 2023
First submitted to journal
12 Sep, 2023

You are reading this latest preprint version

Experimental and Numerical Mudflows Modeling for Runout and Deposition Height Assessmen

Status:

Version 1

Abstract

Figures

1. INTRODUCTION

2. MATHEMATICAL FORMULATION

2.1. Governing equations

2.2. Turbulence model

2.3. Multiphase model

2.3.1. VOF Model

2.3.2. Mixture Model

2.4. Flow viscosity

3. MATERIALS AND METHODS

3.1. Soil characterization

3.1.1. Viscosity as a Newtonian Fluid

3.1.2. Viscosity as a non-Newtonian fluid

3.2. Laboratory channel tests

3.3. Computational dynamic modeling of mudflow

3.3.1. Geometric model and boundary conditions

3.3.2. Mesh generation for the model

3.3.3. Mesh generation for the model

3.4. Case Study

4. RESULTS AND DISCUSSION

4.1. Viscous behavior of the material considering the material as a Newtonian fluid.

4.2. Viscous behavior of the material considering the material as a non-Newtonian fluid.

4.3. Channel tests and numerical models

4.3.1. Numerical model results considering the material as a Newtonian fluid

4.3.2. Numerical model results considering the material as a non-Newtonian fluid

4.3.3. Model Analysis

4.4. Model results in case study: Yangbaodi Landslide

5. CONCLUSIONS AND RECOMMENDATIONS

Declarations

References

Status:

Version 1