Next Article in Journal
Hierarchical Electrode Switching Device Design for Distributed Single-Channel Electrical Resistivity Tomography System
Next Article in Special Issue
Effects of Pile Installation on Existing Tunnels Using Model Test and Numerical Analysis with Medium Density Sand
Previous Article in Journal
Effects of Air Route Alternation and Display Design on an Operator’s Situation Awareness, Task Performance and Mental Workload in Simulated Flight Tasks
Previous Article in Special Issue
Numerical Investigation of Progressive Slope Failure Induced by Sublevel Caving Mining Using the Finite Difference Method and Adaptive Local Remeshing
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Depth-Integrated Two-Phase Modeling of Two Real Cases: A Comparison between r.avaflow and GeoFlow-SPH Codes

1
Civil Engineering For Risk Mitigation, School of Civil, Environmental and Land Management Engineering, Politecnico di Milano, Polo Territoriale di Lecco, Via G. Previati, 1/c, 23900 Lecco, Italy
2
Department of Mathematics and Computer Applied to Civil and Naval Engineering, Universidad Politécnica de Madrid, Calle del Profesor Aranguren, 3, 28040 Madrid, Spain
*
Author to whom correspondence should be addressed.
Submission received: 13 May 2021 / Revised: 11 June 2021 / Accepted: 16 June 2021 / Published: 21 June 2021
(This article belongs to the Special Issue Advanced Numerical Simulations in Geotechnical Engineering)

Abstract

:
Due to the growing populations in areas at high risk of natural disasters, hazard and risk assessments of landslides have attracted significant attention from researchers worldwide. In order to assess potential risks and design possible countermeasures, it is necessary to have a better understanding of this phenomenon and its mechanism. As a result, the prediction of landslide evolution using continuum dynamic modeling implemented in advanced simulation tools is becoming more important. We analyzed a depth-integrated, two-phase model implemented in two different sets of code to stimulate rapid landslides, such as debris flows and rock avalanches. The first set of code, r.avaflow, represents a GIS-based computational framework and employs the NOC-TVD numerical scheme. The second set of code, GeoFlow-SPH, is based on the mesh-free numerical method of smoothed particle hydrodynamics (SPH) with the capability of describing pore pressure’s evolution along the vertical distribution of flowing mass. Two real cases of an Acheron rock avalanche and Sham Tseng San Tsuen debris flow were used with the best fit values of geotechnical parameters obtained in the prior modeling to investigate the capabilities of the sets of code. Comparison of the results evidenced that both sets of code were capable of properly reproducing the run-out distance, deposition thickness, and deposition shape in the benchmark exercises. However, the values of maximum propagation velocities and thickness were considerably different, suggesting that using more than one set of simulation code allows us to predict more accurately the possible scenarios and design more effective countermeasures.

1. Introduction

Every year, hundreds of landslide-related phenomena, including rock avalanches and debris flows, cause loss of human life and the destruction of much else [1]. The risk analysis of this kind of natural hazard requires the evaluation of the stabilization of the source area and the prediction of propagation behavior (runout distance, velocities, thickness, and shape of deposit). In this study, we dealt with the latter issue by using reliable tools for investigating the dynamics of these flows. Once the flow characteristics are predicted, we can develop more precise hazard maps and design more effective countermeasures. There exist many different debris-resisting structures to cope with the threat related to debris avalanches and debris flows, such as bottom drainage screens [2] and artificial barriers [3].
Once unstabilized areas where mass flows are likely to release have been identified, existing programs based on dynamic numerical models are employed to predict flow behaviors and identify possible impact areas. Nowadays, there exists a large number of numerical models [4,5,6] developed to be useful for risk assessment and runout analysis, and the choice of the most suitable simulation code, given the advances in computational capabilities, has raised another important issue. It is always recommended to use more than one simulation to predict specific events or a detailed back-analysis. This allows us to compare numerical results, understand the main differences, and assess potential risks more precisely. Regardless of the selected simulation code’s capabilities, many other factors, including the choice of rheological law [7], the pore–water pressure evolution [8], the presence of erosion [9], and accurate definitions of parameters, affect the results significantly.
This research’s motivation was the significant challenges researchers and users face when attempting to obtain the same results from two different sets of code. We also aimed to evaluate the performance of a depth-integrated, two-phase model implemented in two different sets of code, and investigate the capabilities of these two widely validated sets of code, GeoFlow-SPH [10] and r.avaflow [11], which are based on the same mathematical model have with fundamentally different numerical solvers. The most important thing was to comprehend the methodologies behind the two sets of code, in order to choose an appropriate rheological law and appropriate geotechnical properties, rather than using the same rheological and empirical laws for both sets of code. This is because engineers will use the numerical results obtained for mapping dangerous areas and designing structural countermeasures.
The paper is structured as follows. The first section is dedicated to introducing the two different continuum-based programs and briefly describing the implemented depth-integrated, two-phase mathematical model consisting of the balance equations of mass and linear momentum, and suitable rheological and empirical equations. In the next section, we report the simulations of two different real events. The simulations used the selected sets of code while aiming to evaluate the ability of the depth-integrated two-phase model to reproduce rapid, large-scale landslides. Finally, the obtained numerical results are compared and discussed.

2. Materials and Methods

2.1. Introduction

Nowadays, several advanced simulation tools [4,5,6] are available to define areas in danger of landslides, foresee landslides’ propagation paths, and obtain parameters to design structural and non-structural countermeasures in order to prevent the damage they cause. Here, we have selected GeoFlow-SPH and r.avaflow, which were developed for the propagation modeling of rapid landslides, such as debris flows and rock avalanches, and have been widely validated through the back analysis of real events [11,12,13,14]. GeoFlow-SPH and r.avaflow both employ a two-phase model which uses an interacting solid and fluid mixture; it is similar to the model proposed by Pitman and Le (2005) [15], and based on the depth-integrated mathematical model of Zienkiewicz and Shiomi (1984) [16].
In propagation modeling, the mathematical formulations implemented in sets of code are usually simplified and reduced to two dimensions through depth integration approximations. Since many flow-like landslides can be long and wide compared to their depths, the governing equations can be integrated along the vertical axis. The resulting 2D depth-integrated model provides a good balance between accuracy and time consumption. GeoFlow-SPH and r.avaflow are also based on integrated solutions of the balance equations of linear momentum and mass using the shallow water flow assumption.
We start by explaining how the mentioned simulation programs can be represented. Since the equations are depth-integrated, the particles cannot be defined along the vertical axis. Thus, in depth-integrated models, the flowing mass is divided into a finite number of columns represented by particles. Therefore, two heights corresponding to the solid h s and fluid h w phases will appear in each model by considering the assumption that the fluidized soil is fully saturated, and the sum of these partial heights is the total height h of the mixture. Therefore, the phases are simulated in different layers and interact with each other by interaction force. In Figure 1, the evolution of both phases’ heights with respect to the total height is depicted.
Let n ¯ w = n ¯ and n ¯ s = 1 n ¯ be the fractions of volume V occupied by the water and solid particles, respectively. The heights of the solid phase h s and fluid phase h w are given by:
h s = n ¯ s h h w = n ¯ w h

2.2. r.avaflow Numerical Code

r.avaflow, developed by Mergili et al. (2017, 2018) [11,17], is an open-source computational tool designed for simulating the dynamics of rapid mass flows, including avalanches and debris flows. It represents a GIS-based computational framework, and the sourcing mass and topography can be defined through raster maps. The tool employs the NOC-TVD numerical scheme [18] along with the Pudasaini two-phase model [19], or a simple single-phase, Voellmy-type model. Besides, it considers the interactions between the phases and the entrainment of material along a path. The interested reader will find in the article by Mergili et al. (2017) [11], full details about the r.avaflow numerical code.
In this tool, the flow’s evolution is computed through depth-averaged, two-phase balance equations developed by Pudasaini (2012) [19]. The mass balance equations for the solid and fluid phases, respectively, are given by:
t n s h + x n s h u s + y n s h v s = 0 t n w h + x n w h u w + y n w h v w = 0
The momentum equations for the solid and the fluid phases are written in conservative form, respectively, as:
t n s h u s u s v m + x n s h u s 2 u u s v m + β s x h 2 + y n s h u s v s u v s v m = h S s t n w h u w u w v m + x n w h u w 2 u u w v m + β w x h 2 + y n w h u w v w u v w v m = h S w
where u s v m = γ C u f u s and u f v m = α s / α s C u f u s are the virtual masses for the solid and fluid phases, respectively. γ is the fluid to solid density ratio. C = N v m 1 / n s / n w + γ is the virtual mass force coefficient. N v m = m w / h w m s / h s / u w u s is the virtual mass number. β s = K s g z 1 γ and β w = g z are the hydraulic pressure coefficients for the solid and fluid phases, respectively. K is the earth pressure coefficient and g z the gravity.
The source terms of Equation (3) are as follows:
S s = α s g u s u s tan δ p b s ε p b s b x ε α s γ p b f h x + b x + C D G u f u s u f u s J 1 S w = α w g ε 1 h x h 2 2 p b w + p b w b x 1 α w N R 2 2 u w x 2 + 2 v w y x + 2 u w y 2 X u w ε 2 h 2 + 1 n w N R A 2 x n s x u w u s + y n s x v w v s + n s y u w u s ξ n s u w u s ε 2 n w N R h 2 + 1 γ C D G u f u s u w u s J 1
where g is the non-dimensional gravitational acceleration. δ the internal friction angle. p b s = g z and p b w = p b s / 1 γ are effective basal pressures for the solid and fluid phases. ε is the typical height H to the typical extent of debris flow L . X is the viscous shearing coefficient for fluid and ξ the vertical distribution of the particle concentration. N R = g L ρ w H / n w ρ w and N R A = g L ρ w H / A ρ w are the quasi-Reynolds numbers and the mobility Reynolds numbers, respectively. A is the mobility of fluid at the interface. In the special case in which A = n w , these quasi-Reynolds numbers coincide. C D G is the generalized drag coefficient.
Next, a set of additional functions and laws are needed to complete the conservation of mass and momentum equations presented in the current section. Here, we provide the following items: (1) the basal friction laws, (2) the interaction force laws, and (3) the erosion laws.
Regarding the rheological models implemented in the r.avaflow code, the solid stress is computed based on Mohr-Coulomb plasticity, whereas the fluid stress is modeled as a solid volume-fraction-gradient-enhanced non-Newtonian viscous stress. These models rely on the empirical coefficients C A D for ambient drag and C F F for fluid friction. C A D is a coefficient to be multiplied with the frontal surface and the velocity of the flow to derive air resistance. C F F is Manning’s coefficient used as the fluid friction coefficient, to consider the effects of the roughness of basal surface.
Entrainment is computed through an empirical law in which a user-defined entrainment coefficient E s is multiplied with the total momentum M s + M w at each raster cell. Then, the solid and fluid entrainment rates q E , s and q E , w can be obtained as follows:
q E , s = E s M s + M w n ¯ s , E q E , w = E s M s + M w 1 n ¯ s , E
where n ¯ s , E is the solid fraction of the entrainable material.
In the r.avaflow, the changes of the basal topography are considered by computing the solid and fluid erosion heights, H E , s and H E , w , at each time step, as follows:
H E , s , t = min H E , s t Δ t + q E , s Δ t cos β , H Emax , s H E , w , t = min H E , w t Δ t + q E , w Δ t cos β , H Emax , w
where H Emax , s and H Emax , w are the maximum entrainable depths at the given cell, t is the time passed at the end of the time step, Δ t is the time step length, and β is the local slope of the basal surface.
The solid–fluid interaction is considered by using the two-phase generalized drag model proposed by Pudasaini (2012) [19], a modified Pitman and Le model [15]. It is capable of taking into account the effect of fluid viscosity on drag. The generalized drag coefficient C D G is given by:
C D G = n s n w ρ s ρ w g V T P F R e p + 1 P G R e p m
where V T = g d / γ is the terminal velocity. d the particle diameter. P 0 , 1 is the parameter for a combination of solid and fluid-like contributions to drag resistance. F = γ n w / n s 3 R e p / 180 . R e p = ρ w d u T / η w is the particle Reynolds number. η w the fluid viscosity. G = α f M 1 where M = M R e p varies from 4.65 to 2.4. m is the exponent for drag (1 for linear, 2 for quadratic).
It is important to note that the excess pore pressure terms are disregarded in this model due to the assumption that dilation can be neglected or the permeability is large. Therefore, only hydrostatic fluid pressure is considered: p / z = g z .
Numerical schemes are a fundamental tool for solving differential equations and obtaining approximations of physical problems. In the r.avaflow code, the system of equations is discretized on a staggered grid, in which the system is moved half of the cell size with every time step. Then, the total variation diminishing non-oscillatory central differencing (TVD-NOC) scheme [20,21,22] is employed to solve the model equations. The TVD-NOC has successfully been applied to a large number of mass flow problems [18,21,23,24].

2.3. GeoFlow-SPH Numerical Code

The GeoFlow-SPH code was developed in Madrid by an expert research team for almost a decade. It has previously been applied to theoretical [25], experimental [2], and real case histories [10]. The latest theoretical framework is based on a two-phase, depth-integrated model recently developed by Pastor et al. (2021) [26]. This model is a general approach for two-phase modeling capable of considering many essential physical aspects, such as the phenomena of pore–water pressure evolution and reproducing the dynamic behavior of debris flows by taking into account their soil properties, including permeability and volume stiffness. The two-phase model consists of balance equations of mass and linear momentum expressed, respectively, in quasi-Lagrangian form as follows:
h α t + x α i h α v ¯ α i = n ¯ α e R
ρ α h α d ¯ α v ¯ α d t = ρ α grad P ¯ α 1 2 ρ w b 3 h 2 Δ p ¯ w h grad n ¯ α + τ B α + ρ α b h α + R ¯ h ρ α v ¯ α n ¯ α e R
where sub-index α denotes the solid (s) or fluid (w) phase, n ¯ α is the porosity of the solid or fluid phase ( n ¯ s = 1 n ¯ and n ¯ w = n ¯ ), e R is the erosion rate, and h α = h n ¯ α . b 3 = g is gravitational force and its axis is vertical and points upwards. In the diffusion term of Equation (9), the averaged pressures P ¯ α acting on solid or fluid phases are defined as:
P ¯ s = 1 2 b 3 h h s + Δ p ¯ w h n ¯ ρ s For Solid P ¯ w = 1 2 b 3 h h w Δ p ¯ w h n ¯ ρ w For Fluid
which indicates that, in our special case of interest, total pore–water pressure p w is composed of a hydrostatic part p h y d , which varies linearly from zero at the surface to ρ g h at the bottom, and an excess pore–water pressure Δ p w which should be computed for each time and portion of space.
Next, the mathematical equations given in Equation (9) are completed using a rheological law to compute the basal shear stress τ B . In this study, the numerical analysis was performed through voellmy’s rheological law, which has the same features as the frictional rheological model where the cohesion and viscous terms are disregarded. It is capable of considering the evolution of pore–water pressure at the basal surface. In the case of a pure frictional mass, the basal shear stress is given by:
τ B = ρ d g h Δ p w b v ¯ i v ¯ tan ϕ B + ρ g v ¯ ξ v ¯ i
where h is the propagation height, ϕ B the basal friction angle, v ¯ the depth averaged flow velocity, ξ the turbulence coefficient, and Δ p w b the excess pore–water pressure at the basal surface, which is computed by using a consolidation equation which is given later in this section. It can be seen in the above equation that the basal shear stress τ B will depend on the basal excess pore pressure, and it is modified in accordance with pore–water pressure’s evolution at each node and time step. Take into account that the higher the pore–water pressure, the lower the shear stress at the basal surface.
Regarding entrainment, the empirical erosion law of Hungr has been implemented in the GeoFlow-SPH code. It is based on an algorithm in which the total volume of debris increases in accordance with a specified rate and is given by:
e R = E s h v ¯
which indicates a direct proportionality between the entrainment rate and the product of depth average velocity v ¯ and mobilized soil depth h . E s can be obtained directly from the initial and final volumes of the material and the distance traveled as E s ln V final / V 0 / distance .
The interaction between solid and fluid phases is given by the depth-integrated value of R ¯ . It is provided by the Anderson and Jackson law [27], which is capable of considering large relative velocities. It 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, g is the acceleration of gravity, and m is a constant.
In depth-integrated models, the vertical structure of the magnitudes is lost as the only available information is their depth-integrated values. To overcome this limitation, an additional equation has been implemented in GeoFlow-SPH’s model to describe how pore pressure evolves over time and depth [28]. Here, we recall the consolidation equation describing the evolution of pore pressure along the vertical axis, as follows:
d s Δ p w d t = ρ d b 3 d h d t + c v 2 Δ p w x 3 2 E m 1 1 n ¯ d s n ¯ d t
where ρ d = 1 n ρ s ρ w is the effective density, c v the consolidation coefficient, and E m the odometric modulus. Based on the above equation, the pore–water pressures are also influenced by the propagation height and its velocity, and the porosity variations. The consolidation parameter of c v = k v k w plays an important role, as a high value of the coefficient allows the rapid dissipation of excess pore–water pressure through the consolidation process. In the article by Pastor (2021) [26], the interested reader will find a detailed explanation of the consolidation equation.
Next, we briefly describe a Lagrangian meshless numerical method that is used to discretize the presented depth-integrated equations. Smoothed particle hydrodynamics (SPH) uses it to transform the problems that are basically in the form of partial differential equations (PDEs) to forms suitable for particle-based simulation. SPH was first invented by Lucy (1977) [29] and Gingold and Monaghan (1977) [30] to model astrophysical problems. This technique has been applied in many areas due to its ability to model complicated cases which involve large-displacement deformations, such as the modeling of fast landslides in solid mechanics [4,31,32]. Good reviews can be found in the texts of Liu and Liu (2004) [33] and Li and Liu (2003) [34].
In this model, these 3D problems are transformed into 2D forms by applying a depth-integrated model. As a result, it provides an excellent combination of accuracy and computational time. This technique has been successfully applied to debris flows by Pastor et al. (2014) [10] and Cascini et al. (2014, 2016) [35,36], and debris avalanches by Cuomo et al. (2014, 2016) [37,38].
To discretize the propagating mass in the SPH method, the first step is to present them as a set of nodes, as depicted in the figure, having individual material properties. In this paper, the mathematical equations are discretized while using the two-phase SPH technique to deal with two phases of flow involving solid and water particles. Therefore, two sets of nodes are introduced, one to represent the solid particles’ movement and another to represent the movement of the fluid particles, as depicted in Figure 2.
Then, an interpolation process calculates the relevant properties on each node over neighboring nodes through a kernel function W I J = W x I x J , h without defining any element. The formulation can be expressed as:
f x I = f x J W x I x J , h d x J
where the function f x I is approximated at a position vector x in space.
Once the information is stored into nodes or particles, the integral interpolation of SPH kernel approximation is replaced with a discretized form of summation over all the particles within the region of compact support of kernel function. It is expressed as follows:
f x I h = ˜ N J = 1 f x J W x I x J , h s m ω J
where the infinitesimal volume d x J of the continuous integral representations is replaced by the volume of a neighboring particle ω J .
In the GeoFlow-SPH code, two alternative methods are applied to approximate the pore pressures inside a landslide. In the first alternative, similarly to the most depth-integrated models, the consolidation equation is solved using a quarter-cosine shape function that fulfills boundary conditions with a zero value at the surface and zero gradients at the basal surface to approximate the vertical distribution of pore–water pressure. By neglecting the porosity variations at the basal surface, the following first alternative consolidation equation can be obtained:
d s Δ p w d t = ρ d b 3 d h d t β Δ p w 0 b
where for simplification, we introduce β = c v π 2 / 4 h 2 .
The second alternative numerical model was proposed by Pastor et al. (2015) [13], who combined a finite difference method (FDM) for the 1D equation of a vertical consolidation and depth-integrated SPH model for propagation analysis (SPH-FD model). Pastor et al. (2021) [26] extended this one-phase model to a two-phase model in order to fully approximate the pore pressures inside a landslide. In this technique, a finite-difference mesh, incorporated at each SPH node representing solid particles, is used to discretize pore–water pressures along the vertical axis, as depicted in Figure 3.
The consolidation equation given in Equation (14) can be rewritten in a way that is more suitable for the formulation of FDM, as follows:
d s Δ p w d t = ρ d b 3 d h d t 1 x 3 h + c v 2 Δ p w x 3 2 E m 1 1 n ¯ d s n ¯ d t
where the consolidation equation is discretized by making two changes: (i) height variations and (ii) porosity variations. In this paper, we apply the first alternative model for the cases that are not needed to take into account basal surface permeability and porosity variations of debris flows.

3. Case Studies

This section presents two real case studies and discusses the results obtained from back-analyses performed by the two programs presented in the last sections. The first benchmark exercise was the prehistoric Acheron rock avalanche, which has been back-analyzed by Mergili et al. (2017) [11] using r.avaflow code. This real case was used to assess the validity and performance of the GeoFlow-SPH code. Another validation exercise was the Sham Tseng San Tsuen debris flow, which has been back-analyzed by Pastor et al. (2021) [26] using two-phase SPH modeling implemented in GeoFlow-SPH code. The latter case study was also used to assess the performance of the sets of code by comparing their numerical results.

3.1. Acheron Rock Vvalanche

The Acheron rock avalanche occurred in Canterbury, New Zealand, approximately 1100 years BP [39,40] based on radiocarbon dating analysis. An earthquake may have triggered the landslide, as it is located close to fault zones. In Figure 4, an overview of the Acheron rock avalanche is shown.
The numerical analysis was performed using a 10 m × 10 m digital terrain model. It is estimated that about 6.4 million m 3 of rock mass with the maximum release height of 78.5 m ran down the hillside and traveled 3.5 km with a travel angle of 18° to the south-east. Validation of sets of code relies on the availability of reliable information, including the impact area or deposition area of the event; and estimated deposition thickness, velocities, and final run-out.
The Acheron Rock Avalanche has been back-analyzed by Mergili et al. (2017) [11] using the r.avaflow code. The same back-calculated parameters were used here for both sets of code to simulate the rock avalanche. Flow parameters and coefficients required for both sets of code are given in Table 1.
The densities of solid and fluid particles were taken as ρ s = 2700 kg/m 3 and ρ w = 1000 kg/m 3 , respectively. In the GeoFlow-SPH code, the initial porosity is computed based on densities. On the other hand, the r.avaflow code uses the defined heights of fluid and solid phases to obtain porosities. For both sets of code, an initial porosity of 0.2, for which ρ = 2350 kg/m 3 , was considered.
In the previous numerical analysis [11], the best results were obtained with the dynamic basal friction angle of 17°. In the r.avaflow code, the simulation was based on Mohr–Coulomb-type plastic flow and was performed with an internal angle of friction of 35° and a bed friction angle 17°. On the other hand, the GeoFlow-SPH code uses the voellmy rheological model to analysis case studies. Concerning its rheological parameters, the basal friction angle of 17° and no turbulence coefficient were considered for this particular case. In both sets of code, the fluid phase is modeled based on the Manning formula. However, the fluid friction coefficient was neglected in this study.
The output of r.avaflow is discretized on the basis of geographic information system (GIS) coordinates and consists of raster maps. One of the advantages of the r.avaflow code is performing the simulation model in a GIS in order to manage geo-referenced data and facilitate the use of the results. In the GeoFlow-SPH code, post-processing tools are used to deal with the visualization of output. The main result is the evolution of solid and fluid flow heights, and the outputs are the values of velocity and depositional height computed at each time step. Optionally, total flow pressures and kinetic energies are generated in the r.avaflow code; relative pore–water pressure and total pore–water pressure are generated in the GeoFlow-SPH code.
Figure 5 illustrates the comparison between the numerical results using the r.avaflow and GeoFlow-SPH sets of code, with the parameter values given in Table 1, at different time steps. As shown in the figure, with the same geotechnical parameter values, the two sets of code gave relatively similar deposit distributions, velocities, runout distances, and final deposition thickness values. Both sets of code show the traveling of small portions of the volume in the wrong direction (northern). It should be recalled that the event is an ancient rock avalanche, and perhaps there exist some unrecognizable deposition volumes in the field.
The deposition volume was considered to be the field deposit estimate of 8.9 × 10 6 m 3 for both sets of code using the back-calculated erosion rate given in Table 1. As described in Section 2.2 and Section 2.3, the erosion laws implemented in both sets of code are of empirical type.
During the propagation stage, a significant difference was observed in the numerical results of maximum height evolutions. However, the r.avaflow and GeoFlow-SPH sets of code had similar results regarding the shape of the final deposit and depositional height distribution. The maximum deposit thicknesses calculated by r.avaflow and GeoFlow-SPH sets of code were 31.8 and 27.3 m, respectively.
As shown in Figure 5, both sets of code reached a satisfactory approximation of the depositional area shape. However, the r.avaflow code gave a larger lateral mass spread than GeoFlow-SPH. The excessive lateral spreading might be solved by choosing appropriate combinations of internal friction angle δ and basal friction angel ϕ B . The GeoFlow-SPH code successfully modeled the final deposit, although an isotropic state of stress was assumed, and the earth pressure coefficient (k-value) was considered equal to 1. Lateral spreading also depends on the evolution of pore–water pressure, which the GeoFlow-SPH is capable of considering.
As depicted in Figure 6a, the r.avaflow code overestimated the height values with respect to the GeoFlow-SPH code until 95 s (time of reach). After that, the numerical result of maximum deposit thickness obtained by the r.avaflow decreased during the lateral spreading of the flowing mass until the numerical results of both sets of code became relatively close in terms of the shape of the final deposit and depositional height distribution.
Figure 6b shows the differences in maximum velocity values obtained from the GeoFlow-SPH and r.avaflow programs at the same time steps evaluated during the whole simulation. Both sets of code indicated that maximum velocity peaked at the time step of 20 s and achieved the maximum velocity of about 71 ms 1 . The results show that the r.avaflow and GeoFlow-SPH programs give approximately the same maximum velocities values until the time step of 95 s (time of reach). It is interesting to note that the most significant differences in computed velocities were during the spreading over the deposition area, where the GeoFlow-SPH code predicted that all the particles deposited at the time step of 95 s; on the other hand, at the same time the r.avaflow code indicated that some particles were still traveling at high speed—31 ms 1 .
The numerical analyses clearly show that the obtained calibrated values of rheological parameters can be transferred from one code to another, as both programs are based on the same mathematical formulation, although they employ different empirical laws and numerical solutions.

3.2. Sham Tseng San Tsuen Debris Flow

This debris flow took place in Hong Kong on 23 August 1999. This case was proposed as a benchmarking exercise by the Hong Kong Geotechnical Office with the objective of evaluating the precision of numerical and constitutive models [41]. They also provided a detailed digital map of terrain elevation, including the mobilized mass’s original position and the reservoir’s final position.
The event took place after heavy rain. During the 24 h before the event, the precipitation was 479 mm. The debris flow consisted of 600 m 3 of channeled material, according to the report. It originated from three shallow landslides, of which one was much larger than the others. In fact, the data provided included only this main source. The debris traveled downhill along an incised rocky stream course with negligible debris entrainment, and it destroyed several one-story buildings at the end of the canal. Figure 7 represents a map of the terrain showing the position of landslide A, the track, Sham Tseng San Tsuen village’s location, and the house labeled as “No. 38”, which was destroyed. More information can be found in a geotechnical report provided by Wilson et al. (2005) [42].
This simulation aimed to consider excess pore–water pressures in the two-phase model implemented in GeoFlow-SPH. In the previous simulation, it was assumed that the permeability was very large, and the excess pore pressure terms were neglected.
Densities of soil particles, pore fluid, and mixture were taken as 2400 Kg/m 3 , 1000 Kg/m 3 , and 2000 Kg/m 3 . The drag laws that have been implemented in both sets of code were based on the Anderson and Jackson law [27] model. In this case, the terminal velocity V T and m were taken as 10 4 ms 1 and 1, respectively.
According to the report [42], the triaxial tests show that the material did not present cohesion, and the friction angle varied from 37° to 38°, for which we found tan ϕ B = 0.75 0.78 . Interestingly, the report states that “the distal end of the debris deposition had a travel angle of about 24°”, which is much smaller than the friction angle found in laboratory tests. For the two-phase model developed by Pastor (2021) [26], the values of the friction angle and travel angle coincided well with those obtained in the reports by using the following equation:
tan ϕ a p p = ρ s ρ w ρ s tan ϕ B = 2400 1000 2400 tan 37 ϕ a p p = 24 °
where ϕ B is the friction angle and ϕ a p p the travel angle. The rest of the parameters were similar to those in the previous case study.
Figure 8 presents the propagation predicted by the r.avaflow (right) and GeoFlow-SPH (left) programs at different time steps.
In Figure 8, we can observe the larger velocities of the model simulated in the GeoFlow-SPH code, capable of considering excess pore–water pressures. The main reason for these differences was the pore–water pressure’s evolution in the GeoFlow-SPH code. As described in Section 2.3, the consolidation equation was employed in the two-phase modeling of debris flow to compute the quantity of excess pore–water pressure, which has been implemented in the balance of momentum equation.
Besides, the frictional rheological equation implemented in the numerical code is capable of considering the effects of the pore pressure dissipation and generation. As can be seen in Equation (11), the generation of pore–water pressure is similar in effect to reducing the friction angle, and the dissipation of pore–water pressure is equivalent to increasing the friction angle. To obtain similar numerical results, it is necessary to use a lower value of basal friction angel in the r.avaflow code. Therefore, improvements to the numerical results obtained by the r.avaflow code can be achieved by modifying the values of rheological parameters.
In Figure 9, we provide the numerical results obtained from the GeoFlow-SPH code regarding the basal excess pore–water pressures relative to liquefaction. At the start, the relative pore–water pressure p w r e l was taken as 1, indicating that the flowing material was completely liquefied. Then, the pore–water pressure dissipated during the propagation stage until the debris-flow reached the deposition area.
In the authors’ opinion, the back-analysis of both case studies indicates that the numerical results of r.avaflow are sensitive to variations in the basal friction angle ϕ B and the initial solid fraction n 0 ; and on the other hand, the numerical results of GeoFlow-SPH are more sensitive to the consolidation coefficient c v obtained by using parameters of terminal velocity V T and stiffness of the mixture k v ; consequently, the characteristic times of propagation and consolidation are the key aspects of the GeoFlow-SPH model.
There exist many different structural countermeasures to reduce the impacts of landslides. Bottom drainage screens [2] could be used as an effective control to work against channelized flows with narrow widths, such as the Sham Tseng San Tsuen debris flow case. Artificial barriers and walls [3] could be used as protection against rock and debris avalanches, such as the Acheron rock avalanche, with widespread deposition.

4. Conclusions

This paper described two sets of code, r.avaflow and GeoFlow-SPH, which are similar in mathematical formulation and differ by numerical implementation method. The two programs have been used to back-analyze two real cases that occurred in New Zealand in ancient history and Hong Kong in 1999. In previous works, the former has been back-analyzed with the r.avaflow code to validate the implemented depth-integrated, two-phase model, and the latter has been back-analyzed through the GeoFlow-SPH code to evaluate the two-phase model while considering the space-time evolution of pore–water pressure. The best fit values of geotechnical parameters, which were achieved through one of the programs in the previous works, were used in another program in this study.
The aims of the research were (1) to evaluate the prediction capabilities of programs, which are based on a depth-integrated two-phase model, by comparing the numerical results with available field data, and (2) to investigate whether the same geotechnical parameters calibrated in the one numerical code can be transferred to other numerical sets of code.
The numerical results were compared to the field data, and they evidenced that both sets of code were capable of reproducing the run-out distance, deposition thickness, and shape for the real cases, although some parameters are considered in one numerical code but not in another, such as internal friction angle, earth pressure coefficients, and excess pore pressure evolution. Thus, the results of both sets of code can be considered satisfactory. The first numerical analyses showed that the same values of rheological parameters can be used for different sets of code based on similar mathematical formulae. However, the mathematical formulations used for the second case study were different due to the implementation of a consolidation equation in the Geoflow-SPH code for considering the space-time evolution of pore–water pressure.
There were also considerable differences between the numerical results obtained by the two programs. In the first case study, the GeoFlow-SPH code underestimated the maximum propagation thickness compared to the r.avaflow code. In the second case study, the r.avaflow code underestimated the maximum propagation velocity compared to the GeoFlow-SPH code, which considers the pore-pressure evolution. These results confirm that using more than one simulation code is necessary to predict specific events to assess potential risks and design possible countermeasures.

Author Contributions

S.A.M.T.: conceptualization, formal analysis, writing—original draft, and investigation. S.M.T.: formal analysis, visualization, writing—original draft, and investigation. M.P.: methodology, writing—review and editing, supervision, and resources. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by the Spanish Ministry MINECO under project P-LAND (PID2019-105630GB-I00).

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

All measured data can be obtained by contacting the authors.

Acknowledgments

The authors gratefully acknowledge the support of the Geotechnical Engineering Office, Civil Engineering and Development Department of the Government of the Hong Kong SAR in the provision of the digital terrain models for the Canada landslide case.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Petley, D. Global patterns of loss of life from landslides. Geology 2012, 40, 927–930. [Google Scholar] [CrossRef]
  2. Tayyebi, S.M.; Pastor, M.; Yifru, A.L.; Thakur, V.; Stickle, M.M. Depth Integrated Two-Layer Coupled SPH Models for Debris Flows Simulation, Challenges and Innovations in Geomechanics, IACMAG Challenges Innov. Geomech; Barla, M., Donna, A.D., Sterpi, D., Eds.; Springer International Publishing: New York, NY, USA, 2021; pp. 491–497. ISBN 978-3-030-64518-2. [Google Scholar]
  3. Cuomo, S.; Moretti, S.; D’Amico, A.; Frigo, L.; Aversa, S. Modelling of geosynthetic-reinforced barriers under dynamic impact of debris avalanche. Geosynth. Int. 2020, 1, 65–78. [Google Scholar] [CrossRef]
  4. 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]
  5. Pirulli, M. Numerical Modeling of Landslide Runout. A Continuum Mechanics Approach. Ph.D. Thesis, Politecnico di Torino, Turin, Italy, 2005. [Google Scholar]
  6. MChristen, 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] [Green Version]
  7. Manzanal, D.; Drempetic, V.; Haddad, B.; Pastor, M.; Stickle, M.M.; Mira, P. Application of a New Rheological Model to Rock Avalanches: An SPH Approach. Rock Mech. Rock Eng. 2016, 49, 2353–2372. [Google Scholar] [CrossRef]
  8. 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]
  9. Pirulli, M.; Pastor, M. Numerical study on the entrainment of bed material into rapid landslides. Géotechnique 2012, 62, 959–972. [Google Scholar] [CrossRef]
  10. Pastor, M.; Blanc, T.; Haddad, B.; Petrone, S.; Morles, M.S.; 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, 83, 793–812. [Google Scholar] [CrossRef] [Green Version]
  11. Mergili, M.; Fischer, J.; Krenn, J.; Pudasaini, S.P. r.avaflow v1, an advanced open-source computational framework for the propagation and interaction of two-phase mass flows. Geosci. Model Dev. 2017, 10, 553–569. [Google Scholar] [CrossRef] [Green Version]
  12. 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]
  13. Pastor, M.; Blanc, T.; Haddad, B.; Drempetic, V.; Morles, M.S.; Dutto, P.; Stickle, M.M.; Mira, P.; Merodo, J.A. Depth Averaged Models for Fast Landslide Propagation: Mathematical, Rheological and Numerical Aspects. Arch. Comput. Methods Eng. 2009, 22, 67–104. [Google Scholar] [CrossRef]
  14. 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] [Green Version]
  15. Pitman, E.B.; Le, L. A two-fluid model for avalanche and debris flows. Philos. Trans. R. Soc. A Math. Phys. Eng. Sci. 2005, 20, 1573–1601. [Google Scholar] [CrossRef]
  16. Zienkiewicz, O.C.; Shiomi, T. Dynamic behaviour of saturated porous media: The generalized Biot formulation and its numerical solution. Int. J. Numer. Anal. Methods Geomech. 1984, 8, 71–96. [Google Scholar] [CrossRef]
  17. Mergili, M.; Frank, B.; Fischer, J.; Huggel, C.; Pudasaini, S.P. Computational experiments on the 1962 and 1970 landslide events at Huascarán (Peru) with r.avaflow: Lessons learned for predictive mass flow simulations. Geomorphology 2018, 322, 15–28. [Google Scholar] [CrossRef]
  18. Wang, Y.; Hutter, K.; Pudasaini, S.P. CThe Savage-Hutter theory: A system of partial differential equations for avalanche flows of snow, debris, and mud. Geomorphology 2004, 84, 507–527. [Google Scholar]
  19. Pudasaini, S.P. The Savage-Hutter theory: A general two-phase debris flow model. J. Geophys. Res. Earth Surf. 2012, 117, F03010. [Google Scholar] [CrossRef]
  20. Nessyahu, H.; Tadmor, E. Non-oscillatory central differencing for hyperbolic conservation laws. J. Comput. Phys. 2012, 87, 408–463. [Google Scholar] [CrossRef] [Green Version]
  21. Tai, Y.C.; Noelle, S.; Gray, J.M.; Hutter, K. Shock-Capturing and Front-Tracking Methods for Granular Avalanches. J. Comput. Phys. 2002, 175, 269–301. [Google Scholar] [CrossRef] [Green Version]
  22. Wang, G.; Sassa, K. Factors affecting rainfall-induced flowslides in laboratory flume tests. Géotechnique 2001, 51, 587–599. [Google Scholar] [CrossRef]
  23. Mergili, M.; Schratz, K.; Ostermann, A.; Fellin, W. Physically-based modelling of granular flows with Open Source GIS. Nat. Hazards Earth Syst. Sci. 2012, 12, 187–200. [Google Scholar] [CrossRef] [Green Version]
  24. Pudasaini, S.P.; Krautblatter, M. A two-phase mechanical model for rock-ice avalanches. J. Geophys. Res. Earth Surf. 2012, 119, 2272–2290. [Google Scholar] [CrossRef] [Green Version]
  25. 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]
  26. Pastor, M.; Tayyebi, S.M.; Stickle, M.M.; Yagüe, A.; Molinos, M.; Navas, P.; Manzanal, D. A depth integrated, coupled, two-phase model for debris flow propagation. Acta Geotech. 2021, 1–25. [Google Scholar] [CrossRef]
  27. 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]
  28. Pastor, M.; Quecedo, M.; Merodo, J.A.; Herrores, M.I.; Gonzalez, E.; Mira, P. Modelling tailings dams and mine waste dumps failures. Géotechnique 2002, 52, 579–591. [Google Scholar] [CrossRef]
  29. Lucy, L.B. A numerical approach to the testing of the fission hypothesis. Astron. J. 1977, 82, 1013–1024. [Google Scholar] [CrossRef]
  30. 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]
  31. McDougall, S. A New Continuum Dynamic Model for the Analysis of Extremely Rapid Landslide Motion Across Complex 3D Terrain. Ph.D. Thesis, University of British Columbia, Vancouver, BC, Canada, 2006. [Google Scholar]
  32. 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]
  33. Li, S.; Liu, W.K. Meshfree Particle Methods; Springer: Berlin/Heidelberg, Germany, 2004. [Google Scholar]
  34. Liu, G.R.; Liu, M.B. Smoothed Particle Hydrodynamics: A Meshfree Particle Method; World Scientific: Hackensack, NJ, USA, 2003. [Google Scholar]
  35. Cascini, L.; Cuomo, S.; Pastor, M.; Sorbino, G.; Piciullo, L. SPH run-out modelling of channelised landslides of the flow types. Geomorphology 2014, 214, 502–513. [Google Scholar] [CrossRef]
  36. 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]
  37. Cuomo, S.; Pastor, M.; Cascini, L. and Castorino, G.C. Interplay of rheology and entrainment in debris avalanches: A numerical study. Can. Geotech. J. 2014, 51, 1318–1330. [Google Scholar] [CrossRef]
  38. 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]
  39. Smith, G.M.; Davies, T.R.; McSaveney, M.J.; Bell, D.H. The Acheron rock avalanche, Canterbury, New Zealand—Morphology and dynamics. Landslides 2006, 3, 62–72. [Google Scholar] [CrossRef]
  40. Smith, G.M.; Bell, D.H.; Davies, T.R. The Acheron rock avalanche deposit, Canterbury, New Zealand: Age and implications for dating landslides. N. Z. J. Geol. Geophys. 2012, 55, 375–391. [Google Scholar] [CrossRef]
  41. Pastor, M.; Blanc, T.; Pastor, M.J.; Sánchez, M.; Haddad, B.; Mira, P.; Fernández-Merodo, J.A.; Herreros, I.; Drempetic, V. A SPH depth integrated model with pore pressure coupling for fast landslides and related phenomena. In Proceedings of the 2007 International Forum on Landslide Disaster Management, Hong Kong, China, 10–12 December 2007; Ho, K., Li, V., Eds.; Geotechnical Division, Hong Kong Institution of Engineers: Hong Kong, China, 2017; pp. 987–1014. [Google Scholar]
  42. Wilson, F.M. Report on the Debris Flow at Sham Tseng San Tsuen of 23 August 1999. Findings of the Investigation; GEO Report No. 169; Geotechnical Engineering Office, Civil Engineering and Development Department, The Government of the Hong Kong Special Administrative Region: Hong Kong, China, 2005. [Google Scholar]
Figure 1. Scheme for obtaining partial heights of both phases.
Figure 1. Scheme for obtaining partial heights of both phases.
Applsci 11 05751 g001
Figure 2. SPH interactions for two-phase flow.
Figure 2. SPH interactions for two-phase flow.
Applsci 11 05751 g002
Figure 3. A 1D finite difference mesh at each SPH node representing solid particles.
Figure 3. A 1D finite difference mesh at each SPH node representing solid particles.
Applsci 11 05751 g003
Figure 4. (a) A general view of the source area, (b) the area impacted by the avalanche, and (c) its location (Google Maps).
Figure 4. (a) A general view of the source area, (b) the area impacted by the avalanche, and (c) its location (Google Maps).
Applsci 11 05751 g004
Figure 5. A comparison between two sets of code: the results on the left side were obtained with Geoglow and on the right side with r.avaflow. They show the flow propagation path, with the numerical results of flow-depth and maximum velocity at different times: (a) 0 s, (b) 20 s, (c) 40 s, (d) 60 s, (e) 100 s, and (f) 230 s.
Figure 5. A comparison between two sets of code: the results on the left side were obtained with Geoglow and on the right side with r.avaflow. They show the flow propagation path, with the numerical results of flow-depth and maximum velocity at different times: (a) 0 s, (b) 20 s, (c) 40 s, (d) 60 s, (e) 100 s, and (f) 230 s.
Applsci 11 05751 g005
Figure 6. Comparison of the numerical results of (a) the maximum deposition thickness and (b) the maximum velocity obtained from the r.avaflow and GeoFlow-SPH programs.
Figure 6. Comparison of the numerical results of (a) the maximum deposition thickness and (b) the maximum velocity obtained from the r.avaflow and GeoFlow-SPH programs.
Applsci 11 05751 g006
Figure 7. General view of the landslide location. (Provided by Wilson et al. (2005) [42]).
Figure 7. General view of the landslide location. (Provided by Wilson et al. (2005) [42]).
Applsci 11 05751 g007
Figure 8. A comparison between two sets of code. The results on the left side were obtained with GeoFlow-SPH, and those on the right side with r.avaflow. They show the flow propagation path, with the numerical results of flow-depth and maximum velocity at different times: (a) 0 s, (b) 20 s, (c) 40 s, and (d) 60 s.
Figure 8. A comparison between two sets of code. The results on the left side were obtained with GeoFlow-SPH, and those on the right side with r.avaflow. They show the flow propagation path, with the numerical results of flow-depth and maximum velocity at different times: (a) 0 s, (b) 20 s, (c) 40 s, and (d) 60 s.
Applsci 11 05751 g008
Figure 9. The distribution of relative pore–water pressure at (a) 0 s, (b) 10 s, (c) 20 s, and (d) 60 s.
Figure 9. The distribution of relative pore–water pressure at (a) 0 s, (b) 10 s, (c) 20 s, and (d) 60 s.
Applsci 11 05751 g009
Table 1. Material parameters used in the analysis of the avalanche.
Table 1. Material parameters used in the analysis of the avalanche.
Parameter (Symbol)GeoFlow-SPHr.avaflow
Density of the of solid particle ρ s 2700 kg/m 3 2700 kg/m 3
Density of the of water particle ρ w 1000 kg/m 3 1000 kg/m 3
Mixture density ρ 2350 kg/m 3 -
Internal friction angle δ -35°
Basal friction angle ϕ B 17°17°
Virtual mass force coefficient C -0.5
Terminal velocity V T 0.01 m/s0.01 m/s
Parameter related to drag resistance P -0.5
Particle Reynolds number R e p -1
Exponent for drag -1 D linear, 2 D quadratic- m 11
Quasi-Reynolds numbers N R -30,000
Mobility number N R A -1000
Viscous shearing coefficient for fluid X -0
Solid concentration distribution with depth ξ -0
Ambient drag coefficient C A D -0
Manning’s coefficient C F F 0 s/m 1 / 3 0 s/m 1 / 3
Entrainment coefficient E s 5 × 10 5 m 2 10 7 kg 1
Elastic volumetric stiffness k v 4 × 10 8 N/m 2 -
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Mousavi Tayebi, S.A.; Moussavi Tayyebi, S.; 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. https://0-doi-org.brum.beds.ac.uk/10.3390/app11125751

AMA Style

Mousavi Tayebi SA, Moussavi Tayyebi S, Pastor M. Depth-Integrated Two-Phase Modeling of Two Real Cases: A Comparison between r.avaflow and GeoFlow-SPH Codes. Applied Sciences. 2021; 11(12):5751. https://0-doi-org.brum.beds.ac.uk/10.3390/app11125751

Chicago/Turabian Style

Mousavi Tayebi, Seyed Ali, Saeid Moussavi Tayyebi, and Manuel Pastor. 2021. "Depth-Integrated Two-Phase Modeling of Two Real Cases: A Comparison between r.avaflow and GeoFlow-SPH Codes" Applied Sciences 11, no. 12: 5751. https://0-doi-org.brum.beds.ac.uk/10.3390/app11125751

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