Next Article in Journal
How Does the Spatial Misallocation of Land Resources Affect Urban Industrial Transformation and Upgrading? Evidence from China
Next Article in Special Issue
Inventory and Spatial Distribution of Ancient Landslides in Hualong County, China
Previous Article in Journal
Influence of the 2030 Agenda in the Design of Policies to Fight Poverty and Social Exclusion in Rural and Urban Contexts
Previous Article in Special Issue
A Research on Cohesion Hyperspectral Detection Model of Fine-Grained Sediments in Beichuan Debris Flow, Sichuan Province, China
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Two-Phase Two-Layer Depth-Integrated SPH-FD Model: Application to Lahars and Debris Flows

by
Saeid Moussavi Tayyebi
1,*,
Manuel Pastor
1,
Andrei Hernandez
1,
Lingang Gao
1,2,
Miguel Martin Stickle
1,
Ashenafi Lulseged Yifru
3 and
Vikas Thakur
3
1
Department of Mathematics and Computers Applied to Civil and Naval Engineering, ETS Ingenieros de Caminos, Canales y Puertos, Universidad Politécnica de Madrid, Calle del Profesor Aranguren, 3, 28040 Madrid, Spain
2
College of Water Conservancy and Hydropower Engineering, Hohai University, Nanjing 210098, China
3
Department of Civil and Environmental Engineering, Norwegian University of Science and Technology (NTNU), 7491 Trondheim, Norway
*
Author to whom correspondence should be addressed.
Land 2022, 11(10), 1629; https://doi.org/10.3390/land11101629
Submission received: 24 July 2022 / Revised: 2 September 2022 / Accepted: 14 September 2022 / Published: 22 September 2022

Abstract

:
The complex nature of debris flows suggests that the pore-water pressure evolution and dewatering of a flowing mass caused by the high permeability of soil or terrain could play an essential role in the dynamics behavior of fast landslides. Dewatering causes desaturation, reducing the pore-water pressure and improving the shear strength of liquefied soils. A new approach to landslide propagation modeling considering the dewatering of a mass debris flow has drawn research attention. The problem is characterized by a transition from saturated to unsaturated soil. This paper aims to address this scientific gap. A depth-integrated model was developed to analyze the dewatering of landslides, in which, desaturation plays an important role in the dynamics behavior of the propagation. This study adopted an SPH numerical method to model landslide propagation consisting of pore-water and a soil skeleton in fully or partially saturated soils. In a two-phase model, the soil–water mixture was discretized and represented by two sets of SPH nodes carrying all field variables, such as velocity, displacement, and basal pore-water pressure. The pore-water was described by an additional set of balance equations to take into account its velocity. In the developed two-layer model, an upper desaturated layer and a lower saturated layer were considered to enhance the description of dewatering. This is the so-called two-phase two-layer formulation, which is capable of simulating the entire process of landslides propagation, including the large deformation of soils and corresponding pore-water pressure evolutions, where the effect of the dewatering in saturated soils is also taken into account. A dam-break problem was analyzed through the new and previously developed model. A flume test performed at Trondheim was also used to validate the proposed model by comparing the numerical results with measurements obtained from the experiment. Finally, the model was applied to simulate a real case lahar, which is an appropriate benchmark case used to examine the applicability of the developed model. The simulation results demonstrated that taking into account the effects of dewatering and the vital parameter of relative height is essential for the landslide propagation modeling of a desaturated flowing mass.

1. Introduction

Dewatering occurs when the degree of the saturation of a moving mass is reduced due to the high permeability of soil or terrain that the flow is propagating over. The higher the soil permeability in the propagation mass, the higher the rate of dewatering that can be achieved. It may also play an important in the dynamic behavior of landslides. The dewatering of soil reduces the fluid volume fraction in the upper part of the propagation mass. Consequently, the upper part of the propagation mass may become dry, but the lower part remains fully saturated.
Existing landslide displacement prediction models generally include physical and data-based models [1,2]. The data-based models have recently become popular due to providing a good relation between accuracy and simplicity. However, the physical-based models using continuum dynamics methods are still the only choice to better understand landslide propagation and predict its impacted area.
In some cases, a lahar stops because of basal friction, and its water abandons it, flowing downhill, such as in the case of the Popocatépetl volcano lahar [3,4]. According to the provided record [4], on 22 January 2001, the lahar carrying the pyroclastic material triggered and melted a section of the nearby glacier and flowed down through the drainage gorges. The flow carried 160,000 m 3 of material, including 68 , 000   m 3 of water [5], for approximately 14 km to the village of Santiago Xalitzintla. In this case study, the fluid abandoned the soil skeleton; consequently, the frontal part should have been dominated by the fluid phase. In contrast, the back side of the flow should have been dominated by the dynamic behavior of the solid phase.
The dewatering and desaturation of fast landslides can also be used to reduce the flows’ impact and runout distance [6]. Structure countermeasures, such as bottom drainage screens, can be installed at the flow propagation path to mitigate the potential destructiveness of landslides. In these energy-dissipating structures, negative pore-water pressures were obtained over the surface of the permeable screens [7] or permeable walls [8] due to the dewatering of the flowing mass.
Currently, various continuum dynamics models [9,10,11] exist in advanced simulation tools such as Dan3D [12], RASH3D [13], Ramms [14], r.avaflow [15], and Titan2D [16] to reproduce the propagation stage of fast landslides. In these models, a few realistic simplifying assumptions are considered to obtain the approximate solutions acceptable for many engineering purposes. However, the lack of a unified approach causes several difficulties and uncertainties in an appropriate hazard assessment related to a wide class of phenomena that can occur in both saturated and unsaturated conditions. Therefore, it is essential to develop a generalized computational model to reproduce various types of landslides, taking into account different characteristics, including porosity variations, pore-water pressure evolution, soil permeability, and dewatering effects.
In landslide propagation modeling, it is assumed that, once the flowing mass is deposited, there is no separation of solid and liquid components, as in the case of mudflows. However, it has been reported that some dewatering of the flowing mass occurs following deposition [17,18]. Unlike previous models, this study focused on considering the principles of partially saturated soil mechanics. The crucial part is defining the phenomenon of dewatering that directly affects the dynamic behavior of a flowing mass under the process of the pore-fluid abandoning a solid structure. The developed model is also capable of reproducing the dynamic behavior of landslides propagating over a terrain equipped with structural countermeasures in which dewatering techniques are implemented.
The extended model presented in this work was implemented in the GeoFlow-SPH depth-integrated code created by an expert research team in Madrid almost two decades ago. The implemented model was developed from a one-phase finite element method [19] to a more advanced two-phase smoothed-particle hydrodynamics-finite difference (SPH-FD) method [20] to more properly reproduce the complex dynamics behavior of fast landslides. Several researchers [21,22,23,24] in different institutions have previously applied it to theoretical [4], experimental [25], and real case studies [26].
This paper aims to present some advances in the generalized approach based on the two-phase depth-integrated model proposed by the authors to reproduce landslide propagations with a wide range of soil permeability and densities with a single set of governing equations, and to extend this concept to cases with desaturated soil. The dewatering of a flow and its corresponding change in pore-water pressure can be considered using the proposed model.
Compared with the previous propagation model, the proposed method considering dewatering has many advantages in reproducing landslide propagation with high-permeability soil or terrain. This paper provides a comprehensive review and analysis of the mechanisms and the desaturation implementation techniques.
Following this introduction, we will recall a two-phase model that simulates the fully saturated landslides. This model will be extended to a two-phase two-layer model, which is applicable for landslide cases where the phenomenon of dewatering occurs. The governing equations consist of two linear momentum balance equations, two mass balance equations and a consolidation equation. These equations are completed using a frictional rheological model, an interaction law, and an erosion law. Finally, the model will be applied to three cases where the effects of dewatering should be considered.

2. Material and Methods

2.1. Introduction

First, we need a mathematical model with a set of equations that describes the essential aspects of fast landslide propagation. The depth-integrated two-phase mathematical model, which will be recalled here, is proposed by Pastor et al. [20,27], where the general case of fully saturated soils is considered. It is based on the depth-integrated mathematical model of Zienkiewicz and Mróz [28]. Depth integrated models have been used in hydraulics engineering [29] and frequently applied to landslide propagation [21,30,31,32] due to providing a good relation between accuracy and time consumption.
The model developed by the authors is a two-phase model, similar to those of Pitman and Le [10] and Pudasaini [33], capable of capturing the dynamic of solid and fluid phases. The governing equation consists of two depth-integrated balance equations to describe the propagation phases and a consolidation equation to describe the evolution of excess pore-water pressure along the vertical axis in a depth-integrated model, as follows:
(ⅰ)
Balance equations for mass and linear momentum, respectively:
d ¯ s h d t + . h s v ¯ s = n ¯ s e R For solid phase
d ¯ w h d t + . h w v ¯ w = n ¯ w e R For fluid phase
ρ s h s d s v ¯ s d t = grad 1 2 ρ s b 3 h h s + Δ p w h n ¯ + 1 2 ρ w b 3 h 2 Δ p w h grad n ¯ τ B s + ρ s b h s + h s R ¯ s 1 n ¯ ρ s v ¯ s v ¯ s b e R
ρ w h w d w v ¯ w d t = grad 1 2 ρ w b 3 h h w Δ p w h n ¯ + 1 2 ρ w b 3 h 2 + Δ p w h grad n ¯ τ B w + ρ w b h w + h w R ¯ w n ¯ ρ w v ¯ w v ¯ w b e R
where Δ p w is the excess pore-water pressure. We denote some parameters by the sub-indexes s for solid and w for fluid phases. External forces consist of basal shear stress τ B , gravity force b 3 = g , and interaction force R ¯ . An overbar over a magnitude denotes a depth-integrated value. v ¯ b is the basal slip velocity and e R is an erosion coefficient. As can be seen in Equations (1), (2), (3), and (4), in two-phase models, each phase is characterized by its density ρ , depth-averaged flow velocity v ¯ , height h , and averaged porosity n ¯ . b is the vector of body forces.
(ⅱ)
Consolidation equation:
d s Δ p w d t = ρ ¯ b 3 d h d t + c v 2 Δ p w 0 x 3 2 E m 1 1 n ¯ d s n ¯ d t
which describes pore pressure evolution and consists of three terms, including (i) debris flow height variations, (ii) consolidation along vertical axis, (iii) depth-averaged porosity variations, respectively. c v is a consolidation coefficient obtained by multiplying permeability k w by volumetric stiffness k v . ρ ¯ is the effective density and E m = 3 k v 1 2 ν is the oedometric modulus. ν is the Poisson’s ratio.
This model has previously been applied to the back-analysis of a laboratory-scale flume test [25] performed in Trondheim and several real landslide cases with different types and characteristics, including YuTung channelized debris flow [34], Johnson landscape debris avalanche with widespread deposition [32], Acheron rock avalanche [35], and Sham Tseng San Tsuen debris flow [20].
The porosity variation and dewatering have an important role in saturated soil mechanics. The general depth-integrated model described above allows the porosity to change due to the relative movement between fluid and solid phases. In some cases, a solid flowing material stops because of its high basal friction, but its pore-fluid abandons it. Another interesting example is a debris flow arriving at a permeable retention structure. If the model is applied to these cases, we will obtain a decreasing porosity that can even reach zero, as:
h w 0 and n = h w h s + h w 0
Figure 1a shows a fully saturated material consisting of granular particles and pore-water. A disadvantage of modeling landslides with the depth-averaged equation of motion is their lack of ability to define particles along the depth. To facilitate the interpretation in these models, the flowing mass was divided into a finite number of columns represented by particles, as shown in Figure 1b. Therefore, the total height h of the mixture was divided into two partial heights corresponding to the solid h s and fluid h w phases, which are visualized in two different layers. However, these phases interact with each other through a drag law that depends on the variables of the porosity and permeability of the soil.
Most continuous mediums are composed of various phases. The case of partially saturated soil has three: a mixture consisting of soil grains, water, and air. Additional modifications in propagation modeling are needed to overcome problems involving unsaturated soil and to capture their dynamic behavior and mechanism, which play an important role in the propagation. Therefore, it is essential to consider another material layer for the air phase.
Figure 2a shows a crude simplification of reality, where a saturated layer exists at the bottom, where α is the relative height and α h is the depth of the saturated zone, with high pore-water pressures and a top layer of unsaturated soil. In this model, a degree of saturation S r is considered only for the top layer. Therefore, in this study, we considered a general case where the total height h is divided into three layers, as shown in Figure 2b. The partial heights of solid, fluid, and air can be obtained, respectively, as follows:
h s = 1 n ¯ h
h w = n ¯ α + S r 1 α h
h a = n ¯ 1 S r 1 α h
where n ¯ is the averaged porosity.
Once the porosity reaches a lower bound n ¯ m i n (which depends on the solid grains’ size distribution), the upper part of the sliding mass will be unsaturated. During the computations, we obtained h w and h s from the balance of mass equations; the porosity was computed as:
n ¯ = h w h s + h w n ¯ m i n
If, at a certain time, it reaches the limit, we will compute the depth of the saturated layer as:
α = 1 1 S r 1 n ¯ inf n ¯ inf h w h s S r
where n ¯ n + 1 = h w n + 1 / h n + 1 n ¯ inf .
If we are given h s , h w , S r and α , the total height can be obtained as:
h = h s + h w 1 n ¯ 1 S r 1 α

2.2. Proposed Balance Equations of Mass and Linear Momentum

So far, we have described some particular cases where the soil moisture of the landslide material is lowered from a full-saturation state to a partial-saturation state during the propagation stage. This important phenomenon is called dewatering. Next, the mathematical models and numerical simulations will be developed to take into account both the saturated and unsaturated hydro-mechanical behavior of landslides. The governing equations for both the fully saturated lower layer and the partially saturated upper layer water body were solved using a depth-integrated two-phase model.
Extending the two-phase model proposed by Pastor et al. [20] and considering Equations (1)–(4), the balance equations for solid and fluid phases in both fully and partially saturated soils are written, respectively, as:
(ⅰ)
The balance of mass equations are now written as:
d ¯ s d t 1 n ¯ h + 1 n ¯ h div v ¯ s = 1 n ¯ e R
for the solid phase, and
d ¯ w d t n ¯ α h + n ¯ α h div v ¯ w = n ¯ e R
for the fluid phase.
(ⅱ)
The linear momentum equations are:
ρ s h s d ¯ s v ¯ s d t = grad 1 2 1 n ¯ ρ s h 2 b 3 + n ¯ α h Δ p ¯ w + 1 2 ρ w α 2 h 2 b 3 α h Δ p ¯ w grad n ¯ τ B s + ρ s b h s + h s R ¯ s 1 n ¯ ρ s v ¯ s v ¯ s b e R
for the solid phase, and
ρ w h w d ¯ w v ¯ w d t = grad 1 2 n ¯ ρ w α 2 h 2 b 3 + n ¯ α h Δ p ¯ w + 1 2 ρ w α 2 h 2 b 3 + α h Δ p ¯ w grad n ¯ τ B w + ρ w b h w + h w R ¯ w n ¯ ρ w v ¯ w v ¯ w b e R
for the fluid phase.
The above equations can be written more compactly by introducing:
(ⅰ)
P s and P w defined as:
P s = 1 2 1 n ¯ h 2 b 3 + n ¯ α h Δ p ¯ w ρ s
P w = 1 2 n ¯ α 2 h 2 b 3 + n ¯ α h Δ p ¯ w ρ w
F s and F w defined as:
F s = 1 2 ρ w ρ s α 2 h 2 b 3 α h Δ p ¯ w ρ s
F w = 1 2 α 2 h 2 b 3 + α h Δ p ¯ w ρ w
(ⅱ)
And the source terms:
S s = 1 ρ s h s τ B s + ρ s b h s + h s R ¯ s 1 n ¯ ρ s v ¯ s v ¯ s b e R
S w = 1 ρ w h w τ B w + ρ w b h w + h w R ¯ w n ¯ ρ w v ¯ w v ¯ w b e R
The compacted balance of momentum equations is now written as:
d ¯ s v ¯ s d t = 1 h s grad P s + 1 h s F s grad n ¯ + S s
d ¯ w v ¯ w d t = 1 h w grad P w + 1 h w F w grad n ¯ + S w

2.3. Proposed Consolidation Equation

Pore-water pressure has a significant influence on the initiation and run-out processes. Modeling the evolution in pore-water pressure due to strong coupling between solid and fluid phases is one important part of mechanistic landslide modeling.
In the proposed propagation–consolidation model, the pore-water pressure is allowed to change in space and time from the initiation to the deposition. In this study, numerical analyses were performed in the hypothesis of fully saturated conditions using a depth-integrated two-phase model proposed by Pastor et al. [20], which fits this hypothesis. Then, the analyses were extended to the case of partially saturated conditions caused by dewatering phenomena using the proposed model.
During dewatering, the soil becomes partially saturated, and its dynamic behavior depends on the coupling between the solid skeleton, pore-water, and pore-air. As shown in Figure 3, a two-layer model was proposed. The hydrostatic and excess pore-water pressures were considered for the lower layer. A null pore-water pressure value was assumed at the layer above the water table. In typical cases, air consolidation times are much shorter than those of the water. Therefore, the role of pore-air pressures can be neglected.
The component of the normal stress at the vertical axis σ 3 can be expressed as the sum of the various components, as shown in Figure 3, as follows:
σ 3 = p h y d Δ p w + σ 3
where p h y d = α ρ w b 3 h is the hydrostatic pore pressure with zero air pressure, Δ p w is the excess pore-water pressure, and σ 3 is the effective stress.
The normal stress of layers 1 and 2 can also be obtained as:
σ 3 = α ρ ¯ 1 b 3 h + 1 α ρ ¯ 2 b 3 h
where ρ ¯ 1 and ρ ¯ 2 are the averaged densities of layers 1 and 2:
ρ ¯ 1 = 1 n ρ s + n ρ w
ρ ¯ 2 = 1 n ρ s + n ρ w S r
from which, we obtain the averaged and effective density, respectively, as:
ρ ¯ = α ρ ¯ 1 + 1 α ρ ¯ 2
ρ ¯ = ρ ¯ α ρ w
Finally, the effective stress σ 3 can be obtained as:
σ 3 = b 3 h ρ ¯ d Δ p w
where ρ ¯ d = 1 α ρ ¯ 2 + α ρ ¯ 1 and ρ ¯ 1 = ρ ¯ 1 ρ w .
Suppose we assume that the soil is partially saturated. In that case, when taking into account Equation (5), the excess pore-water evolution can be obtained using the following consolidation equation corresponding to the saturated layer:
d s Δ p w d t = ρ ¯ d b 3 d s h d t + c v 2 Δ p w x 3 2 E m 1 1 n ¯ d s n ¯ d t
The new consolidation equation is combined with the proposed balance equation, an extension of the depth-averaged model proposed by the authors, to predict the motion of a partially saturated debris flow with evolving pore-pressure distributions. The new model enhances the simulation of the dynamics of the pore-water pressure in partially saturated soils.
It is important to note that, in the proposed model, it was assumed that the landslide material consists of solid particles and a pore-water that fully saturates its pores. However, the saturated soil transit to the partially saturated soil during the propagation. In this model, the previous model proposed by the authors was extended to the case of desaturated conditions.
A reliable study of landslide propagation modeling should not ignore these relevant processes, such as pore-water pressure evolution, porosity evolution, and desaturation rate, that affect the dynamics behavior of the flowing mass. Existing models commonly do not consider all of these conditions.

2.4. Rheological and Constitutive Laws

Suitable rheological and constitutive equations must complement the governing equations presented in the previous section. The values of basal shear stress τ B , interaction force R ¯ , and erosion coefficient e R should be computed using these laws.
Dewatering effectively decreases the pore-water pressure in the propagation mass. As a result, the effective stress σ ’ is increased and the regained shear strength reduces the runout of the flow. The mechanisms of dewatering can be explained by the effective stress principle and the friction rheological law. In this study, the solid and fluid phases were modeled using the Voellmy rheological model [36] and the Chezy–Manning formula, respectively, as follows:
τ B s = ρ ¯ d g h Δ p w b v ¯ i v ¯ tan ϕ B + ρ g v ¯ ξ v ¯ i
where h is the total height, ϕ B is the basal friction angle, ξ is the turbulence coefficient, and Δ p w b is the basal excess pore-water pressure obtained using Equation (23), taking into account that, in this rheological model, the basal shear stress depends on the basal excess pore-water pressure, and the higher the pore-water pressure, the lower the shear stress at the basal surface.
τ B w = ρ g n m 2 v ¯ h 7 / 3 v ¯ i
where n m is the Manning’s roughness coefficient.
Two-phase models were applied in both phases: (i) a phase of water and (ii) a phase of solid. These phases interact through a drag law that depends on the porosity and soil permeability variables. There exist different drag laws, such as (i) the Darcy law for the cases with a low relative velocity, (ii) the Anderson–Jackson law [37] for a large relative velocity, and (iii) the Kozeny–Carman law [38,39] for densely packed grains. Pitman and Le [10] used the second law to model the interaction force between the two phases. The interaction force of the law is given by:
R ¯ = n ¯ 1 n ¯ V T n ¯ m ρ s ρ w g v ¯ w v ¯ s
where V T is the terminal velocity and m is related to the Reynold number of the flow (1 for linear and 2 for quadratic). In a recently published paper [40], the authors presented a new solution for the interaction terms of the balance equations, which significantly facilitates the computations of two-phase models.
Regrading the bed entrainment, there exist various empirical [41,42] and process-based [43,44] erosion laws. Pirulli and Pastor [45] briefly described some of these erosion laws. Bed entrainment was not considered in this study. However, the interested reader can find, in the article by Pastor et al. [46], interesting descriptions of an arbitrary Lagrangian–Eulerian (ALE) model implemented in the proposed model to consider pore-water pressure existing in the eroded materials.

2.5. Numerical Models

Concerning numerical models, meshless methods such as the material point method (MPM) [47,48], smooth particle hydrodynamics (SPH) [49,50,51], and other mesh-free methods [52] have been successfully applied in the area of solid mechanics. The obtained partial differential equations should be discretized and transformed to a form suitable for particle-based simulation. In the simulation model developed by Pastor et al. [20], the balance of mass and momentum equations were discretized with the SPH numerical model, whereas the consolidation equation was discretized using a set of finite difference meshes. In this model, two sets of nodes were introduced, one to represent the solid particles’ movement and another to represent the fluid particles’ movement, as well as a 1D finite-difference mesh for describing the pore-water pressure along the vertical axis, as depicted in Figure 4. The meshless numerical method of SPH was invented by Lucy [53] and Gingold and Monaghan [54] to model astrophysical problems. It has also been successfully applied in the area of landslide propagation modeling, particularly in the modeling of debris flows [34,55,56] and other fast landslides [21,31,32] involving large deformations. Concerning depth-integrated SPH models for landslide propagation, it is worth mentioning Pastor et al. [9,57,58,59]. Good reviews can be found in the texts of Li and Liu [60] or Liu and Liu [61].
The results are a set of ordinary differential equations produced in discretized form with respect to time, as follows:
(ⅰ)
Height re-initialization:
h a i = j M a j W i j
where the subscript a denotes the solid (s) or fluid (w) phase. M j is a fictitious volume and W i j is a smoothing kernel function. Subscripts i and j in the above equation denote particles i and j as shown in Figure 4.
(ⅱ)
Linear momentum balance equation:
d ¯ a v ¯ a i d t = N h j = 1 m a j P ¯ a i h a i 2 + P ¯ a j h a j 2 grad W i j N h j = 1 m a j n ¯ a i h a i 2 + n ¯ a j h a j 2 grad W i j + S s
(ⅲ)
Consolidation equation, suitable for the formulation of FDM:
d s Δ p w d t = b 3 d s h d t 1 α ρ ¯ 2 + α ρ ¯ 1 1 x 3 α h + c v 2 Δ p w x 3 2 E m 1 1 n ¯ d s n ¯ d t
The resulting ordinary differential equations are discretized in time with the fourth order Runge–Kutta method [62] in the SPH and the forward time centered space (FTCS) method in the FDM. The presented model is capable of taking into account the evolution of excess pore-water pressure along the saturated layer in a depth-integrated model.

3. Case Studies

This section will apply the coupled two-phase depth-integrated model to describe flows with desaturated soil–water properties. For the sake of simplicity, numerical analyses were performed in the hypothesis of fully saturated conditions in which the saturation degree of flowing material decreases during the propagation stage.

3.1. A Dam Break Problem

One-dimensional dam break exercises provide interesting information and insights into the proposed two-phase two-layer model. In addition, these dam break exercises show the main aspects of the model in simple and controlled situations. The concept of computational models can be comprehended more quickly, showing the results of a simple dam break exercise rather than a landslide over complex terrain. In this study, dam break problems were used to demonstrate the capability of the proposed model to capture the physical aspect of dewatering in fast landslides, which has not been considered before. Figure 5 shows the fully saturated mixture with the partial heights of h s = 6 m and h w = 4 m for the solid and fluid phases, respectively, at the initial conditions. The mixture was assumed to be maintained between two dams, and once the dam on the left side suddenly collapses, the two-phase material flows. The minimum porosity n m i n of the material was assumed to be 0.4. It is important to note that desaturation starts once the porosity decreases and reaches the limit. To consider a fully saturated material, the initial porosity n 0 was assumed to be 0.42, for which, the densities of solid particles, the fluid phase, and bulk were ρ s = 2400 kg / m 3 , ρ w = 1000 kg / m 3 , and ρ w = 1800 kg / m 3 , respectively. The terminal velocity V T was assumed to be equal to 5 m/s, corresponding to a debris flow having high-permeability soil.
Figure 5a,b show the simulations performed through the new two-phase two-layer model considering dewatering effects and the previous model [20], respectively. The numerical results of porosity, partial height, and total height are shown at times 0, 1.5, 3, and 10 s. As shown in Figure 5a1,b1, the initial porosity n 0 of the flowing mass is 0.42. It increases at the frontal part of the flow and decreases at the backside during the propagating stage. Once the porosity reaches its minimum, the degree of saturation starts decreasing, and, consequently, the phenomenon of dewatering, which the previous model [20] was not able to reproduce, occurs at a rate that depends on the soil permeability. The total height and porosity profiles exhibit a substantial evolution due to the high soil permeability and desaturation rate. As depicted in Figure 5a4, the porosity reaches unity at the frontal part of the flow once the water leaves the solid skeleton. The solid particles cannot be mobilized at this stage due to low interaction forces between phases. Therefore, the material behavior can be described as drained due to the rapid dissipation of pore-water pressures and the short consolidation time. This saturated granular flow can be representative of wet rock avalanches or debris flows consisting of soils with a high permeability.

3.2. Flume Test

Flow velocity can be reduced using structure countermeasures. Some of these structures, such as flexible nets and bottom drainage screens, are designed to dissipate basal pore-water pressure and desaturate the flowing mass. In return, soil particles regain contact friction, so the shearing resistance of the moving mass increases.
In this simulation, the capacity of the extended model was evaluated to reproduce a debris flow on a small-scale flume test. The laboratory test was performed by Yifru [63] at Trondheim to illustrate the effectiveness of the structural countermeasure of the bottom drainage screen. For the first time, the mentioned flume test was simulated using the 2D depth-integrated SPH-FD model. The small-scale laboratory test was performed in a 10 m long channel. It has two parts: (i) the main channel, which is 6 m long and 0.3 m wide with a 17 inclination, and (ii) a deposition area, 4 m long and 2.2 m wide with a 2 inclination. Figure 6a shows the initial configuration of the test with fully saturated debris material with a total volume of 25 L and a solid concentration of 60%. Considering that the source material’s solid concentration is 60%, the densities of soil particles, interstitial fluid, and mixture were taken with values of 2750 kg / m 3 , 1000 kg / m 3 , and 2050 kg / m 3 , respectively. The interested reader will find a detailed description of the experiments in the text by Yifru [63].
Figure 6b shows that the debris flow is passing through the boundary wall at the main channel. As shown in Figure 6c,d, once the flowing material crosses the bottom drainage screen, the velocity of the solid particles starts declining rapidly. At this stage, the basal total pore-water pressure dissipates during the propagation over the screen. Consequently, the excess pore-water pressure instantly becomes equal, but opposite in sign, to the hydrostatic pressure at the basal surface. Due to the low permeability of the soil, only a considerable amount of water at the basal surface abandons the solid skeleton. Therefore, the pore-water pressure still exists in the body of the flow, and the presented depth-integrated two-phase SPH-FD is capable of considering this essential physical aspect. It is important to note that the change in concentration of the debris flow due to the drainage of pore-water through the basal screen causes the soil particles to regain their contact friction. Therefore, the basal shear stress of the flowing material increases, which the presented model is capable of considering when using the Voellmy Equation (33).
Figure 7 shows the numerical results of the final deposit thickness consisting of three layers of solid, fluid, and air at a propagation time of 10 s produced by the proposed model. As shown in the figure, there are two distinguishing parts of the debris deposit. (i) The first is the debris mass accumulated over the permeable screen with a deposit thickness of 7 cm. (ii) The second is the debris mass accumulated with a maximum height of 2 cm in the deposition area. The numerical results satisfactorily reproduce the runout distance and final deposition heights of the laboratory test over the permeable screen and the deposition area. The desaturation rate is relatively low in this case due to a high water content and low soil permeability. However, a significant amount of material trapped over the screen shows the effectiveness of this energy-dissipating structure against fast landslides.
The flow heights recorded by sensors (fh-1, fh-2, and fh-4) and the computed flow heights (SPH-1, SPH-2, and SPH-4) at distances of 2.6 m, 5 m, and 5.6 m are plotted in Figure 8. The arrival time of the flow front obtained from the presented numerical model fits with the laboratory results of the fh-1 and fh-2. The flow shows an exception at the fh-4, where the simulation looks faster than what was observed in the laboratory. This could be because of the complexities of a depth-integrated model used to reproduce a small-scale flume test equipped with boundary walls. In such cases, the velocity distribution along the channel width has a maximum at the center, and zero at the sides, which was not considered in the presented model. However, the depth-averaged flow height was able to reproduce the flow heights and deposition thicknesses of the flume test well.
Out of curiosity, we assumed that the permeability of the soil used in the laboratory was high. Thus, we simulated it with the presented method. Figure 9 shows the deposition shape and thickness from two perspectives obtained from the numerical model. Such a flume test containing high-permeability soil has not been physically modeled. As we expected in this scenario, almost the whole material is trapped over the screen. The reason is that, the higher the soil permeability in the propagation mass, the higher the desaturation rate that can be achieved. Consequently, the lower the pore-water pressure, the higher the shear strength at the base of the material.
The primary purposes of installing such structures are to separate some amount of water from the flowing mass, weaken the pore-water pressure, and decrease the propagating mass’s velocity. To capture the realistic behavior of solid and fluid phases and their interaction, the models must take into account pore-water pressure evolution and desaturation to simulate flow-like landslides on a terrain equipped with permeable walls or screens.

3.3. Real Case: Popocatépetl Volcano Lahar

In some cases, the phenomena of flowing mass dewatering may occur following deposition [17,18], such as the Popocatépetl volcano lahar in which the flowing mass stopped because of basal friction, and its water abandoned it, flowing downhill [3,4]. In this case, the lahar behaved as a debris flow at the initial stages of the propagation. Once the phenomenon of dewatering occurred, it behaved as a debris flood throughout its propagation up to the deposition.
Popocatépetl is an active volcano with destructive eruptive phases, located 70 km southeast of Mexico City with a population of 24 million. Hundreds of lahar events have occurred in the last two decades along the Huiloac channel. These lahars took place in 1997, 2001, and 2017 due to the earthquakes or the melting of the Ventorillo glacier, which triggered two large lahar and hundreds of shallow landslides. The volume of water in the mixture that caused the 2001 lahar flow was much greater than the volume corresponding to the 1997 and 2017 debris flows.
It has been reported [4] that:
“On 22 January 2001 the destruction of a dome inside the crater triggered a pyroclastic flow that crossed the glacier and flowed down through the drainage gorges. The flow abraded the glacier and triggered the melting and release of 717 , 000 m 3 of water, based on calculations from photogrammetry records of the reduction in the mass of the glacier. Four hours later a lahar was initiated in the valley head areas; this carried the pyroclastic flow materials down Huiloac Gorge for approximately 14 k m to the immediate surroundings of Santiago Xalitzintla. The lahar transported 160,000 m 3 of material, including 68,000 m 3 of water [5]. This lahar behaved as a debris flow throughout its course [64]. Initial lahar volume was estimated by different methods [5,65] leading to oscillating values from 1.85 to 3.30 × 10 5 m 3 for 1997 lahar and 1.57 to 2.44 × 10 5 m 3 for 2001 lahar."
The densities of the soil particles, pore-fluid, and mixture were taken with values of 2400 kg / m 3 , 1000 kg / m 3 , and 2000 kg / m 3 , respectively. These values agree with the case description of the flow mass consisting of a considerate portion of water (i.e., over 68,000/160,000 = 42.5% of the flow mass is water). Entrainment was not considered in the simulation due to the observation [4] that this lahar event had a low erosive power. Due to the high permeability of the material [66], the terminal velocity V T and oedometric modulus E m were taken as 15 × 10 3 m / s and 10 7 N / m 2 , respectively.
The Popocatépetl lahar was simulated by Haddad et al. [4] using a one-phase depth-integrated SPH model. They assumed the material to be of Bingham type, its yield strength τ y , and viscosity μ being selected as 540 Pa and 15 Pa.s, respectively. In one-phase models, the Bingham model could be a better choice than a frictional rheological model to simulate a slurry debris flood. To apply a one-phase frictional rheological model to the lahar with the watery nature, an extremely low friction angle should be adopted to reproduce the lahar’s runout distance, which does not agree with the reported observations. In the one-phase model, the results are more sensitive to the parameter of the frictional angle. On the other hand, employing more input parameters in a two-phase model allows for better matching between the simulation results and measurements. In addition, in the frictional rheological model, the basal shear stress depends on the basal pore-water pressure, and it is modified following pore-water pressure evolution at each node and time step. The higher the pore-water pressure, the lower the shear stress at the basal surface. Therefore, the apparent friction angle introduced would be much smaller than the measured one because vertical consolidation could have changed it during the propagation phase. In this simulation, the solid and fluid phases were modeled using the frictional rheological model and the Chezy–Manning formula, respectively. Concerning its rheological parameters, the basal friction angle of 11 , the Manning roughness coefficient of 0.015 s / m 1 / 3 , and no turbulence coefficient were considered.
The simulation of the Popocatétl lahar was performed using the proposed two-phase two-layer model for the entire propagation of the event. Depth-integrated numerical tools for simulating this event are ideal due to its large runout zone and debris volume. Figure 10 summarizes the simulation results and shows that the proposed model successfully simulated the lahar.
Figure 10a shows the topography and the initial condition used in the model to simulate the Popocatépetl lahars. The numerical analysis was performed using a 15 m × 15 m digital terrain model. The initial mass was schematized into two sets of 12,275 SPH computational points, one representing solid particles and the other representing fluid particles, 5 m spaced. Increasing the number of computational elements could improve numerical results, with the trade-off being an increase in computational time, particularly in this case study with a large domain area.
We present in Figure 10b the isolines of the basal Péclet number ( P e ) at time 25 s, where h 1 m , ρ w = 1000 Kg m 3 , g 10 m s 2 and E m 10 7 Pa . This results in P e 10 2 , which indicates that the consolidation time is much larger than the dewatering time.
As shown in Figure 10c, a large volume of fluid has abandoned the solid skeleton at time 500 s. Figure 10d,e show the propagation heights of the solid and fluid phases, respectively, at the time 1500 s, indicating that the fluid phase at the front of the flow is the dominant phase. Figure 10f shows that the simulated lahar reaches the final deposition, and that the computed and measured run-out distance is very close to each other (14 km). Until now, no numerical simulation model has been reported that is capable of capturing the mechanisms of dewatering in an actual landslide during the propagation stage.
It is important to note that the numerical results of the developed model are more sensitive to the variables of the terminal velocity, stiffness of the mixture, and, consequently, the consolidation coefficient. Therefore, the characteristic times of propagation and consolidation are the key aspects of the proposed model.

4. Discussion

In the previous section, we performed three simulations using the proposed depth-integrated two-phase two-layer SPH-FD model. The simulation results of the flow heights, porosity profiles, basal Péclet number, and runout distance were shown at different time steps. In all of the case studies, the porosities of the mixtures evolve during the propagating stage. The degree of saturation remains equal to unity until the porosity decreases and reaches the minimum porosity of the material. After this point, desaturation starts at a rate that depends on the permeability of the soil. In these cases, a large volume of fluid abandons the solid skeleton due to the high permeability of soil or terrain. The results show that the fluid dominates the frontal part of the flowing mass with a high soil permeability whereas the solid dominates the backside. In such cases, the interaction forces between the phases are not strong enough to mobilize the solid particles. For these reasons, dewatering has been used as a means to mitigate the propagation of fast landslides. Therefore, applying the proposed model that takes into account dewatering is essential to capture the realistic behavior of solid and fluid phases and their interaction. The proposed two-layer model considering dewatering in a flowing mass can be implemented in all depth-integrated models regardless of the numerical method considered.
These numerical tools have many limitations, particularly depth-integrated models due to disregarding vertical accelerations. We can only discuss here the capabilities of these computational models in capturing some physical aspects of landslide propagation. Consequently, engineers or users must verify that these tools with many restrictions are capable of carrying out their case studies.

5. Conclusions

The purpose of this paper is to present a simple depth-integrated model that describes the propagation of fast landslides containing desaturated materials caused by the high permeability of soil or terrain. It is a two-phase model capable of discretizing a soil–water mixture in two sets of SPH nodes carrying all field variables, such as the velocity, displacement, and basal pore-water pressure. It is also a two-layer model, with an upper desaturated layer and a lower saturated layer, which enhances the description of dewatering. Pore-water pressure evolution can be described along the vertical axis in a depth-integrated model using a consolidation equation discretized by the finite difference method.
The model can be applied to two cases where the phenomenon of dewatering may occur. Dewatering can be the consequence of (i) pore-fluid leaving the solid skeleton due to the high permeability of the soil or (ii) pore-fluid draining through a structural countermeasure, such as a bottom drainage screen.
First, a one-dimensional dam break exercise was used to demonstrate the main feature of the proposed model compared to the previous model. In this exercise, the limitations of the previous model were shown. Then, the model was validated by reproducing a complex dynamic behavior of a debris flow on a flume test equipped with a bottom drainage screen by comparing the numerical results with the measurements obtained from the experiment. Finally, the model was applied to a real lahar case where the phenomena of flowing mass dewatering occurred.
These cases demonstrated that it is essential to take into account (i) the effects of desaturation by using a two-layer model and (ii) the effects of dewatering in propagation modeling by implementing the vital parameter of relative height in the governing equations. In fact, by comparing the numerical results obtained for landslide problem modeling to the experimental-measured and field-measured results, we improved the capability of previous models to simulate the desaturated soil behavior for landslide propagation analysis.
The proposed two-phase two-layer SPH-FD model is capable of correctly performing the time–space evolution of pore-water pressures while considering desaturation during the whole propagation stage over an impermeable and permeable bottom boundary. Using this model, it is possible to reproduce the effect of desaturation on the drained behavior of saturated coarse-grained soil.

Author Contributions

Conceptualization, S.M.T. and M.P.; methodology, S.M.T. and M.P.; software, S.M.T. and M.P.; validation, S.M.T.; formal analysis, S.M.T.; investigation, A.H. and L.G.; resources, A.L.Y. and V.T.; data curation, S.M.T.; writing—original draft preparation, S.M.T.; writing—review and editing, M.P.; visualization, S.M.T.; supervision, M.P. and M.M.S.; project administration, M.P. and M.M.S.; funding acquisition, M.P. and M.M.S. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by the Inter-American Development Bank (Project number: ES-T1343) and the Ministerio de Ciencia e Innovación (Grant Number: PID2019-105630GBI00).

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

The data that support the findings of this study are available on request from the corresponding author. The data are not publicly available.

Acknowledgments

The authors gratefully acknowledge the support of the Department of Continuum Mechanics and Theory of Structures, ETSI Caminos, Canales y Puertos, Universidad Politecnica de Madrid, and especially Angel Yague, Diego Manzanal, Miguel Molinos, and Pedro Navas. Authors would like to thank the administrative and technical support of the ETSI Caminos, Canales y Puertos at the Universidad Politécnica de Madrid.

Conflicts of Interest

The authors declare no conflict of interest.

Abbreviations

The following notations are used in this manuscript:
b 3 Gravity force b 3 = g [m/s 2 ]
b Vector of body forces[m/s 2 ]
c v Consolidation coefficient[m 2 /s]
E m Oedometric modulus[N/m 2 ]
e R Erosion coefficient[m/s]
hPropagation height[m]
k v Elastic volumetric stiffness[N/m 2 ]
k w Permeability tensor[m 4 /Ns]
MFictitious volume[m 3 ]
mExponent for drag (1 for linear and 2 for quadratic)[-]
m a Mass. The subscript a denotes the solid (s) or fluid (w) phase.[m 3 ]
n ¯ Averaged porosity[-]
n m Manning roughness coefficient[s/m 1 / 3 ]
R ¯ Interaction force[N/m 3 ]
S r Degree of saturation[%]
V T Terminal velocity[m/s]
v ¯ Averaged velocity[m/s]
v ¯ b Basal slip velocity[m/s]
WSmoothing kernel function[-]
α Relative height[-]
Δ p w Excess pore-water pressure[Pa]
Δ p w b Excess pore-water pressure at the basal surface[Pa]
ν Poisson’s ratio[-]
ξ Turbulence coefficient[m/s 2 ]
ρ Bulk density[kg/m 3 ]
ρ ¯ Averaged density[kg/m 3 ]
ρ s Density of the of solid particle[kg/m 3 ]
ρ w Density of the of water particle[kg/m 3 ]
ρ ¯ Effective density ρ ¯ = ρ α ρ w [kg/m 3 ]
σ 3 Normal stress tensor at the vertical axis[N/m 2 ]
σ 3 Effective stress tensor at the vertical axis[N/m 2 ]
τ B Basal shear stress[N/m 2 ]
ϕ B    Basal friction angle[ ]

References

  1. Ma, J.; Xia, D.; Guo, H.; Wang, Y.; Niu, X.; Liu, Z.; Jiang, S. Metaheuristic-based support vector regression for landslide displacement prediction: A comparative study. Landslides 2022, 19, 2489–2511. [Google Scholar] [CrossRef]
  2. Ma, J.; Xia, D.; Wang, Y.; Niu, X.; Jiang, S.; Liu, Z.; Guo, H. A comprehensive comparison among metaheuristics (MHs) for geohazard modeling using machine learning: Insights from a case study of landslide displacement prediction. Eng. Appl. Artif. Intell. 2022, 114, 105150. [Google Scholar] [CrossRef]
  3. Coviello, V.; Capra, L.; Norini, G.; Dávila, N.; Ferrés, D.; Márquez-Ramírez, V.H.; Pico, E. Earthquake-induced debris flows at Popocatépetl Volcano, Mexico. Earth Surf. Dyn. 2021, 9, 393–412. [Google Scholar] [CrossRef]
  4. Haddad, B.; Pastor, M.; Palacios, D.; Muñoz-Salinas, E. A SPH depth integrated model for Popocatépetl 2001 lahar (Mexico): Sensitivity analysis and runout simulation. Eng. Geol. 2010, 114, 312–329. [Google Scholar] [CrossRef]
  5. Muñoz-Salinas, E.; Castillo-Rodríguez, M.; Manea, V.; Manea, M.; Palacios, D. Lahar flow simulations using LAHARZ program: Application for the Popocatépetl volcano, Mexico. J. Volcanol. Geotherm. Res. 2009, 182, 13–22. [Google Scholar] [CrossRef]
  6. Yifru, A.L.; Laache, E.; Norem, H.; Nordal, S.; Thakur, V.K. Laboratory investigation of performance of a screen type debris-flow countermeasure. HKIE Trans. 2018, 25, 129–144. [Google Scholar] [CrossRef]
  7. Gonda, Y. Function of a debris-flow brake. Int. J. Eros. Control Eng. 2009, 2, 15–21. [Google Scholar] [CrossRef]
  8. Li, X.; Zhao, J.; Kwan, J.S. Assessing debris flow impact on flexible ring net barrier: A coupled CFD-DEM study. Comput. Geotech. 2020, 128, 103850. [Google Scholar] [CrossRef]
  9. McDougall, S.; Hungr, O. A model for the analysis of rapid landslide motion across three-dimensional terrain. Can. Geotech. J. 2004, 41, 1084–1097. [Google Scholar] [CrossRef]
  10. Pitman, E.B.; Le, L. A two-fluid model for avalanche and debris flows. Philos. Trans. R. Soc. A Math. Phys. Eng. Sci. 2005, 363, 1573–1601. [Google Scholar] [CrossRef]
  11. Pudasaini, S.P.; Hutter, K. Avalanche Dynamics: Dynamics of Rapid Flows of Dense Granular Avalanches; Springer: Berlin/Heidelberg, Germany, 2007. [Google Scholar] [CrossRef]
  12. Hungr, O.; McDougall, S. Two numerical models for landslide dynamic analysis. Comput. Geosci. 2009, 35, 978–992. [Google Scholar] [CrossRef]
  13. Pirulli, M. Numerical Modeling of Landslide Runout. A Continuum Mechanics Approach. Ph.D. Thesis, Politecnico di Torino, Torino, Italy, 2005. [Google Scholar]
  14. Christen, M.; Kowalski, J.; Bartelt, P. RAMMS: Numerical simulation of dense snow avalanches in three-dimensional terrain. Cold Reg. Sci. Technol. 2010, 63, 1–14. [Google Scholar] [CrossRef]
  15. Mergili, M.; Jaboyedoff, M.; Pullarello, J.; Pudasaini, S.P. Back calculation of the 2017 Piz Cengalo–Bondo landslide cascade with r.avaflow: What we can do and what we can learn. Nat. Hazards Earth Syst. Sci. 2020, 20, 505–520. [Google Scholar] [CrossRef]
  16. Sheridan, M.; Stinton, A.; Patra, A.; Pitman, E.; Bauer, A.; Nichita, C. Evaluating Titan2D mass-flow model using the 1963 Little Tahoma Peak avalanches, Mount Rainier, Washington. J. Volcanol. Geotherm. Res. 2005, 139, 89–102. [Google Scholar] [CrossRef]
  17. Pierson, T.C. Erosion and deposition by debris flows at Mt Thomas, North Canterbury, New Zealand. Earth Surf. Process. 1980, 5, 227–247. [Google Scholar] [CrossRef]
  18. Rodine, J.D.; Johnson, A.M. The ability of debris, heavily freighted with coarse clastic materials, to flow on gentle slopes. Sedimentology 1976, 23, 213–234. [Google Scholar] [CrossRef]
  19. Pastor, M.; Quecedo, M.; Gonzalez, E.; Herreros, M.I.; Merodo, J.A.; Mira, P. Modelling of Landslides: (II) Propagation. In Degradations and Instabilities in Geomaterials; Darve, F., Vardoulakis, I., Eds.; Springer: Vienna, Austria, 2004; pp. 319–367. [Google Scholar] [CrossRef]
  20. Pastor, M.; Tayyebi, S.M.; Stickle, M.M.; Yagüe, Á.; Molinos, M.; Navas, P.; Manzanal, D. A depth integrated, coupled, two-phase model for debris flow propagation. Acta Geotech. 2021, 16, 2409–2433. [Google Scholar] [CrossRef]
  21. Cascini, L.; Cuomo, S.; Pastor, M.; Rendina, I. SPH-FDM propagation and pore water pressure modelling for debris flows in flume tests. Eng. Geol. 2016, 213, 74–83. [Google Scholar] [CrossRef]
  22. Lin, C.; Pastor, M.; Yague, A.; Tayyebi, S.M.; Stickle, M.M.; Manzanal, D.; Li, T.; Liu, X. A depth-integrated SPH model for debris floods: Application to Lo Wai (Hong Kong) debris flood of August 2005. Géotechnique 2019, 69, 1035–1055. [Google Scholar] [CrossRef]
  23. Krušić, J.; Abolmasov, B.; Marjanović, M.; Pastor, M.J.; Tayyebi, S. Numerical modeling of Selanac debris flow propagation using SPH code. In Proceedings of the 13th International Symposium on Landslides, Cartagena, Colombia, 22–26 February 2021. [Google Scholar]
  24. Longo, A.; Pastor, M.; Sanavia, L.; Manzanal, D.; Martin Stickle, M.; Lin, C.; Yague, A.; Tayyebi, S.M. A depth average SPH model including μ(I) rheology and crushing for rock avalanches. Int. J. Numer. Anal. Methods Geomech. 2019, 43, 833–857. [Google Scholar] [CrossRef]
  25. Tayyebi, S.M.; Pastor, M.; Yifru, A.L.; Thakur, V.K.S.; Stickle, M.M. Two-phase SPH–FD depth-integrated model for debris flows: Application to basal grid brakes. Géotechnique 2021, 1–16. [Google Scholar] [CrossRef]
  26. Pastor, M.; Blanc, T.; Haddad, B.; Petrone, S.; Sanchez Morles, M.; Drempetic, V.; Issler, D.; Crosta, G.B.; Cascini, L.; Sorbino, G.; et al. Application of a SPH depth-integrated model to landslide run-out analysis. Landslides 2014, 11, 793–812. [Google Scholar] [CrossRef]
  27. Pastor, M.; Yague, A.; Stickle, M.M.; Manzanal, D.; Mira, P. A two-phase SPH model for debris flow propagation. Int. J. Numer. Anal. Methods Geomech. 2018, 42, 418–448. [Google Scholar] [CrossRef]
  28. Zienkiewicz, O.C.; Mróz, Z. Generalized Plasticity Formulation and applications to Geomechanics. In Mechanics of Engineering Materials; Pearson: London, UK, 1984; pp. 655–679. [Google Scholar]
  29. Saint-Venant, A. Theorie du mouvement non permanent des eaux, avec application aux crues des rivieres et a l’introduction de marees dans leurs lits. Comptes-Rendus l’Académie Des. Sci. 1871, 73, 147–273. [Google Scholar]
  30. Savage, S.B.; Hutter, K. The dynamics of avalanches of granular materials from initiation to runout. Part I: Analysis. Acta Mech. 1991, 86, 201–223. [Google Scholar] [CrossRef]
  31. Cuomo, S.; Pastor, M.; Capobianco, V.; Cascini, L. Modelling the space–time evolution of bed entrainment for flow-like landslides. Eng. Geol. 2016, 212, 10–20. [Google Scholar] [CrossRef]
  32. Tayyebi, S.M.; Pastor, M.; Stickle, M.M.; Yagüe, Á.; Manzanal, D.; Molinos, M.; Navas, P. Two-phase SPH modelling of a real debris avalanche and analysis of its impact on bottom drainage screens. Landslides 2022, 19, 421–435. [Google Scholar] [CrossRef]
  33. Pudasaini, S.P. A general two-phase debris flow model. J. Geophys. Res. Earth Surf. 2012, 117, F03010. [Google Scholar] [CrossRef]
  34. Tayyebi, S.M.; Pastor, M.; Stickle, M. Two-phase SPH numerical study of pore-water pressure effect on debris flows mobility: Yu Tung debris flow. Comput. Geotech. 2021, 132, 103973. [Google Scholar] [CrossRef]
  35. Tayebi, S.A.M.; Tayyebi, S.M.; Pastor, M. Depth-Integrated Two-Phase Modeling of Two Real Cases: A Comparison between r.avaflow and GeoFlow-SPH Codes. Appl. Sci. 2021, 11, 5751. [Google Scholar] [CrossRef]
  36. Voellmy, A. Uber die Zerstorungskraft von Lawinen. Schweiz. Bauztg. 1955, 73, 159–162. [Google Scholar]
  37. Anderson, T.B.; Jackson, R. Fluid Mechanical Description of Fluidized Beds. Equations of Motion. Ind. Eng. Chem. Fundam. 1967, 6, 527–539. [Google Scholar] [CrossRef]
  38. Kozeny, J. Ueber kapillare Leitung des Wassers im Boden. Sitzungsber Akad. Wiss. Wien 1927, 136, 271–306. [Google Scholar]
  39. Carman, P.C. Flow of Gases through Porous Media; Butterworths: London, UK, 1956. [Google Scholar]
  40. Tayyebi, S.M.; Pastor, M.; Stickle, M.M.; Yagüe, Á.; Manzanal, D.; Molinos, M.; Navas, P. SPH numerical modelling of landslide movements as coupled two-phase flows with a new solution for the interaction term. Eur. J. Mech. B/Fluids 2022, 96, 1–14. [Google Scholar] [CrossRef]
  41. Hungr, O.; McDougall, S.; Bovis, M. Entrainment of material by debris flows. In Debris-Flow Hazards and Related Phenomena; Springer: Berlin/Heidelberg, Germany, 2005; pp. 135–158. [Google Scholar] [CrossRef]
  42. Egashira, S. Mechanism of Sediment Erosion and Deposition of Debris Flow. Jpn. Soc. Eros. Control Eng. 1993, 46, 45–49. [Google Scholar]
  43. Iverson, R.M. Elementary theory of bed-sediment entrainment by debris flows and avalanches. J. Geophys. Res. Earth Surf. 2012, 117, F03006. [Google Scholar] [CrossRef]
  44. Vicari, H.; Tran, Q.A.; Nordal, S.; Thakur, V. MPM modelling of debris flow entrainment and interaction with an upstream flexible barrier. Landslides 2022, 19, 2101–2115. [Google Scholar] [CrossRef]
  45. Pirulli, M.; Pastor, M. Numerical study on the entrainment of bed material into rapid landslides. Géotechnique 2012, 62, 959–972. [Google Scholar] [CrossRef]
  46. Pastor, M.; Tayyebi, S.M.; Stickle, M.M.; Molinos, M.; Yague, A.; Manzanal, D.; Navas, P. An Arbitrary Lagrangian Eulerian (ALE) finite difference (FD)-SPH depth integrated model for pore pressure evolution on landslides over erodible terrains. Int. J. Numer. Anal. Methods Geomech. 2022, 46, 1127–1153. [Google Scholar] [CrossRef]
  47. Molinos, M.; Martín Stickle, M.; Navas, P.; Yagüe, Á.; Manzanal, D.; Pastor, M. Toward a local maximum-entropy material point method at finite strain within a B-free approach. Int. J. Numer. Methods Eng. 2021, 122, 5594–5625. [Google Scholar] [CrossRef]
  48. Cuomo, S.; Perna, A.D.; Martinelli, M. Material point method (MPM) hydro-mechanical modelling of flows impacting rigid walls. Can. Geotech. J. 2021, 58, 1730–1743. [Google Scholar] [CrossRef]
  49. Lin, C.; Wang, X.; Pastor, M.; Zhang, T.; Li, T.; Lin, C.; Su, Y.; Li, Y.; Weng, K. Application of a Hybrid SPH - Boussinesq model to predict the lifecycle of landslide-generated waves. Ocean Eng. 2021, 223, 108658. [Google Scholar] [CrossRef]
  50. Chalk, C.; Pastor, M.; Peakall, J.; Borman, D.; Sleigh, P.; Murphy, W.; Fuentes, R. Stress-Particle Smoothed Particle Hydrodynamics: An application to the failure and post-failure behaviour of slopes. Comput. Methods Appl. Mech. Eng. 2020, 366, 113034. [Google Scholar] [CrossRef]
  51. Lian, Y.; Bui, H.H.; Nguyen, G.D.; Zhao, S.; Haque, A. A computationally efficient SPH framework for unsaturated soils and its application to predicting the entire rainfall-induced slope failure process. Géotechnique 2022, 1–19. [Google Scholar] [CrossRef]
  52. Navas, P.; Molinos, M.; Stickle, M.M.; Manzanal, D.; Yagüe, A.; Pastor, M. Explicit meshfree u -pw solution of the dynamic Biot formulation at large strain. Comput. Part. Mech. 2022, 9, 655–671. [Google Scholar] [CrossRef]
  53. Lucy, L.B. A numerical approach to the testing of the fission hypothesis. Astron. J. 1977, 82, 1013–1024. [Google Scholar] [CrossRef]
  54. Gingold, R.A.; Monaghan, J.J. Smoothed particle hydrodynamics—Theory and application to non-spherical stars. Mon. Not. R. Astron. Soc. 1977, 181, 375–389. [Google Scholar] [CrossRef]
  55. Cascini, L.; Cuomo, S.; Pastor, M.; Sorbino, G.; Piciullo, L. SPH run-out modelling of channelised landslides of the flow type. Geomorphology 2014, 214, 502–513. [Google Scholar] [CrossRef]
  56. Cuomo, S.; Pastor, M.; Cascini, L.; Castorino, G.C. Interplay of rheology and entrainment in debris avalanches: A numerical study. Can. Geotech. J. 2014, 51, 1318–1330. [Google Scholar] [CrossRef]
  57. Pastor, M.; Haddad, B.; Sorbino, G.; Cuomo, S.; Drempetic, V. A depth-integrated, coupled SPH model for flow-like landslides and related phenomena. Int. J. Numer. Anal. Methods Geomech. 2009, 33, 143–172. [Google Scholar] [CrossRef]
  58. Rodriguez-Paz, M.; Bonet, J. A corrected smooth particle hydrodynamics formulation of the shallow-water equations. Comput. Struct. 2005, 83, 1396–1410. [Google Scholar] [CrossRef]
  59. Nikooei, M.; Choi, C.E. Towards Depth-Averaged Modelling of the Decay of Granular Flows by Deposition. Comput. Geotech. 2022, 148, 104792. [Google Scholar] [CrossRef]
  60. Li, S.; Liu, W.K. Meshfree Particle Methods, 1st ed.; Springer: Berlin/Heidelberg, Germany, 2004. [Google Scholar] [CrossRef]
  61. Liu, G.R.; Liu, M.B. Smoothed Particle Hydrodynamics: A Meshfree Particle Method; World Scientific: Singapore, 2003; p. 472. [Google Scholar] [CrossRef]
  62. Hirsch, C. Numerical Computation of Internal & External Flows: Fundamentals of Numerical Discretization; John Wiley & Sons, Inc.: New York, NY, USA, 1988. [Google Scholar]
  63. Yifru, A.L. Investigation of a Screen Structure for Mitigating Debris-Flows Along Coastal Roads. Ph.D. Thesis, Norwegian University of Science and Technology, Trondheim, Norway, 2020. [Google Scholar]
  64. Capra, L.; Poblete, M.; Alvarado, R. The 1997 and 2001 lahars of Popocatépetl volcano (Central Mexico): Textural and sedimentological constraints on their origin and hazards. J. Volcanol. Geotherm. Res. 2004, 131, 351–369. [Google Scholar] [CrossRef]
  65. Sheridan, M.F.; Hubbard, B.; Bursik, M.I.; Abrams, M.; Siebe, C.; Macias, J.L.M.; Delgado, H. Gauging short-term volcanic hazards at Popocatépetl. Eos Trans. Am. Geophys. Union 2001, 82, 185–189. [Google Scholar] [CrossRef]
  66. Haddad, B.; Palacios, D.; Pastor, M.; Zamorano, J.J. Smoothed particle hydrodynamic modeling of volcanic debris flows: Application to Huiloac Gorge lahars (Popocatépetl volcano, Mexico). J. Volcanol. Geotherm. Res. 2016, 324, 73–87. [Google Scholar] [CrossRef]
Figure 1. Scheme for obtaining partial heights of both phases. h s and h w are the partial heights corresponding to the solid and fluid phases.
Figure 1. Scheme for obtaining partial heights of both phases. h s and h w are the partial heights corresponding to the solid and fluid phases.
Land 11 01629 g001
Figure 2. Magnitudes characterizing the partially saturated debris flow components.
Figure 2. Magnitudes characterizing the partially saturated debris flow components.
Land 11 01629 g002
Figure 3. Pore-water pressure distribution in solid–fluid mixture of total height h. The lower layer (Layer 1) is fully saturated, whereas the upper layer (Layer 2) could be partially saturated or dry S r = 0 .
Figure 3. Pore-water pressure distribution in solid–fluid mixture of total height h. The lower layer (Layer 1) is fully saturated, whereas the upper layer (Layer 2) could be partially saturated or dry S r = 0 .
Land 11 01629 g003
Figure 4. SPH integration in the particle support for two-phase: (1) soil–soil (I-J) and (2) soil–water (I-K) with a 1D finite-difference mesh at each SPH node that represents a solid particle.
Figure 4. SPH integration in the particle support for two-phase: (1) soil–soil (I-J) and (2) soil–water (I-K) with a 1D finite-difference mesh at each SPH node that represents a solid particle.
Land 11 01629 g004
Figure 5. The profiles of the propagation heights and their corresponding porosity of dam break simulations performed through (a) the new and (b) the previous model [20], at times 0, 1.5, 3, and 10 s.
Figure 5. The profiles of the propagation heights and their corresponding porosity of dam break simulations performed through (a) the new and (b) the previous model [20], at times 0, 1.5, 3, and 10 s.
Land 11 01629 g005
Figure 6. Three-dimensional depth-integrated modeling of a flume test equipped by a bottom drainage screen at times (a) 0 s, (b) 1.0 s, (c) 1.5 s and (d) 6.0 s.
Figure 6. Three-dimensional depth-integrated modeling of a flume test equipped by a bottom drainage screen at times (a) 0 s, (b) 1.0 s, (c) 1.5 s and (d) 6.0 s.
Land 11 01629 g006
Figure 7. The numerical results of final deposit thickness consisting of solid, fluid, and air phases.
Figure 7. The numerical results of final deposit thickness consisting of solid, fluid, and air phases.
Land 11 01629 g007
Figure 8. Experimental results vs. numerical results of the flow heights.
Figure 8. Experimental results vs. numerical results of the flow heights.
Land 11 01629 g008
Figure 9. The numerical results of deposition thickness of flume test, containing high permeability soil, in two different views: (a) isometric Z and (b) plane XZ.
Figure 9. The numerical results of deposition thickness of flume test, containing high permeability soil, in two different views: (a) isometric Z and (b) plane XZ.
Land 11 01629 g009
Figure 10. (a) The terrain and the initial mass, (b) basal Péclet number at time 25 s, (c) overtaking the liquid phase over the solid phase, (d) the propagation height of the solid phase (m) at the time 1500 s, (e) the propagation height of the fluid phase (m) at the time 1500 s, and (f) the lahar reaching the final deposition with the computed run-out distance (m).
Figure 10. (a) The terrain and the initial mass, (b) basal Péclet number at time 25 s, (c) overtaking the liquid phase over the solid phase, (d) the propagation height of the solid phase (m) at the time 1500 s, (e) the propagation height of the fluid phase (m) at the time 1500 s, and (f) the lahar reaching the final deposition with the computed run-out distance (m).
Land 11 01629 g010
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Tayyebi, S.M.; Pastor, M.; Hernandez, A.; Gao, L.; Stickle, M.M.; Yifru, A.L.; Thakur, V. Two-Phase Two-Layer Depth-Integrated SPH-FD Model: Application to Lahars and Debris Flows. Land 2022, 11, 1629. https://doi.org/10.3390/land11101629

AMA Style

Tayyebi SM, Pastor M, Hernandez A, Gao L, Stickle MM, Yifru AL, Thakur V. Two-Phase Two-Layer Depth-Integrated SPH-FD Model: Application to Lahars and Debris Flows. Land. 2022; 11(10):1629. https://doi.org/10.3390/land11101629

Chicago/Turabian Style

Tayyebi, Saeid Moussavi, Manuel Pastor, Andrei Hernandez, Lingang Gao, Miguel Martin Stickle, Ashenafi Lulseged Yifru, and Vikas Thakur. 2022. "Two-Phase Two-Layer Depth-Integrated SPH-FD Model: Application to Lahars and Debris Flows" Land 11, no. 10: 1629. https://doi.org/10.3390/land11101629

APA Style

Tayyebi, S. M., Pastor, M., Hernandez, A., Gao, L., Stickle, M. M., Yifru, A. L., & Thakur, V. (2022). Two-Phase Two-Layer Depth-Integrated SPH-FD Model: Application to Lahars and Debris Flows. Land, 11(10), 1629. https://doi.org/10.3390/land11101629

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop