**4. Conclusions**

Good corporate governance is primarily the responsibility of every company and regulations at the national level taking into consideration specifics of the national economy, and the latest developments of governance practices and regulations at the European or even global level should ensure that certain governance standards are respected. Therefore, it is important that both hard law and soft law (i.e. especially corporate governance codes) provide comprehensive corporate governance framework, thereby encouraging the introduction of high governance standards and best practices in the companies' corporate governance system. This is of key importance for performance, growth and long-term sustainability of companies.

The findings of research studies and analyses discussed in this contribution show that the Slovenian CG Code has been playing an important role in developing corporate governance practice in Slovenia. Especially this is true for Slovenian-listed companies supporting cognitions by Aguilera and Cuervo-Cazurra [31] about the issuer's ability to enforce changes in the companies' governance system. Codes that are developed by the stock markets have strongest enforceability, since they are designed as the norm of operation, and thus having a greater impact on the promotion of good governance. The CG Code itself as well as mandatory disclosure of compliance with the code's recommendations serves as a guideline to different groups of stakeholders by clearly describing the particularities of the Slovene business world. Disclosures in the CG Statement based on 'comply or explain' approach should be specific and of high quality so that shareholders, investors and other stakeholders get a transparent and a reliable picture of the company's governance system.

The LJSE analyses from 2012 and 2015 of disclosures of compliances with the Slovenian CG Code [40, 50] show that the number of specific, high-quality explanations of deviations describing besides deviations also alternative solutions has increased. Even though these results indicate that companies have become aware of the importance of good disclosure practices, their share is still relatively low and therefore improvements are needed in this respect. Such situation is not specific for Slovenia, but has been noticed in other European countries where companies often use standard explanations of deviations, see [10]. The 'comply or explain' approach is effective only when high-quality explanations of deviations are disclosed. That is a way we find of crucial importance to raise awareness of the companies' key stakeholders on the main features of the high-quality explanations. According to the EC Recommendations, the high-quality explanations of deviations mean [32, Article 18]

• avoiding the use of the standardized language,

A closer look at the results shows that the share of companies disclosing deviations from particular recommendations has increased in the observed period. Half of the statistically most

Another research on corporate governance codes in Slovenia was conducted as a part of research on measuring the corporate governance quality by applying *the SEECGAN Index of Corporate Governance* [38, 39]. The main goal of the SEECGAN Index of Corporate Governance (the SEECGAN Index) is to enable assessment of the quality of governance practice in individual companies as well as the situation regarding overall corporate governance practice in a national economy. Additionally, it makes possible to compare corporate governance practices in South Eastern European countries since the SEECGAN Index methodology pays attention to the particularities of these economies. Seven governance categories are assessed [51, 52]:

All companies of the prime and standard market, in total 22 companies, that were listed on the Ljubljana Stock Exchange in June 2014 were explored. The main source of data was annual reports for the year 2013. Additionally, reports and documents published on the companies' websites were analysed. Research results revealed that more than three-quarter of the prime and the standard market companies disclosed a code. The majority of companies referred to one of the official codes. All standard market companies and 88.9% of the prime market companies disclosed the compliance with the corporate governance code and explained the deviations from it. Even though the disclosure of compliance with the chosen code is obligatory for the prime and the standard market companies in Slovenia [45, 46], one of the prime

Good corporate governance is primarily the responsibility of every company and regulations at the national level taking into consideration specifics of the national economy, and the latest developments of governance practices and regulations at the European or even global level should ensure that certain governance standards are respected. Therefore, it is important that both hard law and soft law (i.e. especially corporate governance codes) provide comprehensive corporate governance framework, thereby encouraging the introduction of high governance standards and best practices in the companies' corporate governance system. This is of

key importance for performance, growth and long-term sustainability of companies.

market companies did not disclose compliance with the code [39].

frequent deviations were those from the recommendations on the transparency.

**1.** Structure and Governance of Boards,

76 Corporate Governance and Strategic Decision Making

**4.** Corporate Social Responsibility,

**6.** Corporate Risk Management and

**7.** Compensation/Remuneration.

**5.** Audit and Internal Control,

**3.** Shareholders' Rights,

**4. Conclusions**

**2.** Transparency and Disclosure of Information,

EC recommends companies [32; Section III, paragraph 8, 33; Section III, paragraph 8] to clearly state which specific recommendations they do not comply with and for each deviation, they should

The research findings show that the level of compliance with the codes' recommendations has been increasing in Slovenia. However, as stated in both reports [40, 50], we cannot make a firm conclusion on the actual level of compliance with the CG Code's recommendations. Companies may not disclose all deviations or may find them as unimportant in their attempt to disclose compliance with as many recommendations as possible. That is a way companies should be made aware of the main purpose of the corporate governance code and accompanying 'comply or explain' approach since 'departing from a provision in the code could in some cases allow a company to govern itself more effectively' [32, p. 44]. A non-compliance with the 'best practice' which is accompanying with an explanation of how the alternative approach achieves the goal of the non-adopted recommendation can present significant benefit when creating the governance system that best suits the company's specific circumstances, see [36]. Companies should be aware of the flexibility enabled by the 'comply or explain' approach, and develop a governance system that in the best possible way addresses the company's specifics. The practice of not disclosing all deviations could be a very dangerous one since it can raise doubt about the implementations of the rest recommendations for which a company disclose compliance, see [40].

Both analyses on disclosures of compliances with the CG Code [40, 50] provide important cognitions on the adoption of the CG Code in Slovenian companies. Findings of such analyses reveal improvements in the governance practice and indicate those areas where changes are required. That is a way such monitoring and analysis should be done on the regular basis. Since we can observe high concentrated ownership in Slovenia [50] and companies with controlling shareholders are found to be less prone to disclose information [26], a regular monitoring of disclosures is of great importance. The EC recommends that public or specialized bodies should regularly monitor corporate governance statements published by companies in order to make 'comply or explain' approach effective [32, Article 19]. Shareholders should also perform effective monitoring in order to encourage good-quality explanations [32, Article 20]. Shareholders should play an active role in improving the quality of explanations in Slovenia as well. A dialogue between shareholders, a management board and a supervisory board is of great importance in the process of creating a suitable governance system. External institutions as professionals in monitoring the quality of disclosures [40] cannot do this. However, such professionals can play an important role in the process of monitoring due to knowledge and expertise they possess.

Reporting on the monitoring results can considerably contribute to better understanding of the code's recommendations among companies, promote debate and thus foster awareness of the underlying issues, see [26]. Regular monitoring of the codes adoption can provide legislators, policy makers and stock exchanges with an important insight into the effectiveness and efficiency of the codes, thus providing basis for developing and updating the recommendations 'in order to address the potential failures of corporate governance mechanisms' [10, p. 222]. Such monitoring can be the opportunity for regulators to 'make the rules less ambiguous' [26, p. 196] as it is the case with the last revision of the Slovenian CG Code that considered the findings of analysis of disclosures of compliances with the Slovenian CG Code from 2009 for the 2011–2014 periods.

The research studies and analysis not only in Slovenia but also in other surroundings deal mainly with the disclosures of compliance with the codes' recommendation. However, we believe that future research should address not only the statements on (non)-compliances but also how companies implement and practice the code's recommendations. Effective governance is demonstrated by the implementation of the regulations and recommendation in the practice and 'whether a code really contributes to improving practices depends on the extent to which companies actually comply with the recommendations in the code and to what extent compliances leads to changes in corporate behaviour' [12, p. 63]. We believe that a more in-depth analysis of the declared and implemented governance arrangements and their consequences is needed. An important contribution in this direction can be the SEECGAN Index that enables to measure how the codes' recommendations and national regulations contribute to the quality of corporate governance practice. The SEECGAN Index enables the comparisons of governance practices among South Eastern European countries, thereby creating a platform for debate about the best governance practices considering the specifics on national economies in that part of Europe.

Future research should also address the relationship among the code's compliance and the company's performance in general as well as in Slovenia. None of the researches and analyses conducted in Slovenia have addressed this question yet. We find this issue to be of great importance especially since mixed results about the impact of the level of compliance with the code's recommendations on companies' performance can be found in the literature, see [10]. Some research studies even showed that financial performance could justify non-compliance [7]. Diverse and mixed results can be explained by cognitions of several authors that corporate governance is a complex construct influenced by many factors, see [37]. Both the research and the practice regulated by different forms of hard and soft laws should adequately address the complexity of corporate governance construct. We hope that findings presented in this paper contribute to better understanding on how the codes of good governance as a form of soft law address this complexity and where improvements are required.


**ADVANCEMENTS AND BREAKTHROUGHS IN ULTRASOUND IMAGING**

Edited by **Gunti Gunarathne**

#### **Advancements and Breakthroughs in Ultrasound Imaging**

http://dx.doi.org/10.5772/46053 Edited by Gunti Gunarathne

#### **Contributors**

Philip Njemanze, G P P Gunarathne, Prakash S. Hiremath, Prema T. Akkasaligar, Sharan Badiger, Fenster, Chih-Yen Tu, Wei-Chih Liao, Chuen-Ming Shih, Wu-Huei Hsu, Hung-Jen Chen, Chia-Hung Chen, Frank Lindseth, David Romero-Laorden, Timm Westhoff, Nikolaos Pagonas, Stergios Vlatsas, Steven L. Orebaugh

#### **© The Editor(s) and the Author(s) 2013**

The moral rights of the and the author(s) have been asserted.

All rights to the book as a whole are reserved by INTECH. The book as a whole (compilation) cannot be reproduced, distributed or used for commercial or non-commercial purposes without INTECH's written permission. Enquiries concerning the use of the book should be directed to INTECH rights and permissions department (permissions@intechopen.com).

Violations are liable to prosecution under the governing Copyright Law.

Individual chapters of this publication are distributed under the terms of the Creative Commons Attribution 3.0 Unported License which permits commercial use, distribution and reproduction of the individual chapters, provided the original author(s) and source publication are appropriately acknowledged. If so indicated, certain images may not be included under the Creative Commons license. In such cases users will need to obtain permission from the license holder to reproduce the material. More details and guidelines concerning content reuse and adaptation can be foundat http://www.intechopen.com/copyright-policy.html.

#### **Notice**

Statements and opinions expressed in the chapters are these of the individual contributors and not necessarily those of the editors or publisher. No responsibility is accepted for the accuracy of information contained in the published chapters. The publisher assumes no responsibility for any damage or injury to persons or property arising out of the use of any materials, instructions, methods or ideas contained in the book.

First published in Croatia, 2013 by INTECH d.o.o. eBook (PDF) Published by IN TECH d.o.o. Place and year of publication of eBook (PDF): Rijeka, 2019. IntechOpen is the global imprint of IN TECH d.o.o. Printed in Croatia

Legal deposit, Croatia: National and University Library in Zagreb

Additional hard and PDF copies can be obtained from<EMAIL_ADDRESS>
Advancements and Breakthroughs in Ultrasound Imaging Edited by Gunti Gunarathne p. cm. ISBN 978-953-51-1159-7 eBook (PDF) ISBN 978-953-51-7159-1


A. Dumitrache

128 Nonlinearity, Bifurcation and Chaos – Theory and Applications

onset of flow bifurcations [40, 41].

**6. Conclusion**

neglected \* *y* 1 .

distribution of shear stress.

empirical constant was need.

examined.

The nonlinear flow equations, the initial solution used for the numeric computations and the length of the airfoil midpart with a small or zero curvature are the principal factors for the

In the beginning of the chapter we have achieved an analytic solution that approximates a two-dimensional Coanda flow. The validity of the results is limited to the cases given by *b R*/ 1, since in the tangential component of momentum equation the curvature was

The validity of the laminar solution is for Reynolds numbers smaller than the critical value which is of the order 4 3 10 . For the turbulent flow the similarity assumption has been reduced primarily to the general development along the flow, leaving the transverse

In this model the edge of the intermittent flow has been neglected, so that assuming constant turbulent viscosity in cross section has lead to similar velocity profiles for both the laminar and turbulent flow. However, the general configuration of the flow development was different. In order to specify the development of this self-modeled flow only one

In many applications that use boundary layer control by tangential blowing, the solid surface downstream of blowing slot is strongly curved and, in this case, the prediction of jet involves separation and a more accurate knowledge of the flow (radial and tangential

The compressible Reynolds-averaged Navier–Stokes (RANS) equations have been solved for circulation control (CC) airfoil flows. Different turbulence models have been considered for closure, including the Spalart–Allmaras model with and without a curvature correction and the shear stress transport (SST) model of Menter. Numerical solutions have been computed with a structured grid solver. The effect of mesh density on the solutions has been

We have investigated the characteristics of various Coanda surfaces, involving smooth curved surfaces and a polygonal curved surface with flap. Using the FLUENT code we have

Further, we have taken interest in the detailed behavior of an existing Coanda ejector model, used in propulsion systems. For numerical investigations we have used an implicit formulation of RANS equations for axisymmetric flow with a shear stress transport *k*

(SST model) turbulence model. The numerical results have been obtained for a total pressure range of 1-5 bars, imposed at the reservoir inlet. The goal was to investigate the influence of various geometric parameters and pressure ratios on the Coanda ejector performance. The effect of various factors, such as the pressure ratio, primary nozzle and ejector configurations on the system performance has been evaluated based on the

analyzed the distribution pressure and separation on the considered surfaces.

pressure - velocity profiles) which has been done by CFD methods.

*Institute of Statistics and Applied Mathematics of the Romanian Academy, Bucharest, Romania*


**8. Role of mentee**

**The driver of relationship**: He identifies the skills, knowledge and/or goals they want to achieve and communicates them to their mentor.

**Development planner**: Mentee works with his/her mentor on deficiencies and seeks resources for learning by identifying people and information.

**Contributor**: He looks for opportunities to give back to their mentor.

**Life-long learner**: The standards of the Liaison Committee on Medical Education (LCME) require the educational program must include instructional opportunities for active learning and independent study to foster the skills necessary for lifelong learning. By enrolling in the mentoring program, the mentee will gain the traits of life-long learner [20].


[11] Glowinski R, Pan T-W, Hesla TI, Joseph DD. A distributed lagrange multiplier/fictitious domain for particulate flows. The International Journal of Multiphase Flow. 1999;**25**: 755-794. DOI: 10.1016/S0301-9322(98) 00048-2

[12] Peskin CS. The immersed boundary method. Acta Numerica. 2002;**11**:479-517. DOI: 10.1017/S0962492902000077

[13] Yu Z, Shao X. A direct-forcing fictitious domain method for particulate flows. Journal of Computational Physics. 2007;**227**:292-314. DOI: 10.1016/j. jcp.2007.07.027

[14] Wachs A. A DEM-DLM/FD method

for direct numerical simulation of particulate flows: Sedimentation of polygonal isometric particles in a Newtonian fluid with collisions. Computers & Fluids. 2010;**38**(8): 1608-1628. DOI: 10.1016/j. compfluid.2009.01.005

[15] Davis RH. In: Tory EM, editor. Velocities of Sedimenting Particles in Suspension, in: Sedimentation of Small Particles in Viscous Fluid, Advances in Fluid Mechanics. Southampton: Computational Mechanics Publications; 1996. pp. 161-198

[16] Guillou S, Makhloufi R. Effect of a shear-thickening rheological behaviour on the friction coefficient in a plane channel flow: A study by direct numerical simulation. Journal for Non-Newtonian Fluid Mechanics. 2007;**144**: 73-86. DOI: 10.1016/J.JNNFM.2007. 03.008

[17] Verjus R, Guillou SS, Ezersky A, Angilella JR. Chaotic sedimentation of particle pairs in a vertical channel at low Reynolds number: Multiple states and routes to chaos. Physics of Fluids. 2016;**28**:123303. DOI: 10.1063/ 1.4968559

[18] Verjus R. Contribution à la modélisation des processus de sédimentation-consolidation dans les ports et les estuaires: Etude bi-échelle. [thesis]. University of Nancy; 2015

[19] Yu Z, Phan-Thien N, Fan Y, Tanner RI. Viscoelastic mobility problem of a system of particles. Journal for Non-Newtonian Fluid Mechanics. 2002;**104**: 87-124. DOI: 10.1016/S0377-0257(02) 00014-9

[20] Happel J, Brenner H. Low Reynolds Number Hydrodynamics, 552P. The Hague: Martinus Nijhof; 1973

[21] Jayaweera KOLF, Mason BJ. The behaviour of freely falling cylinders and cones in a viscous fluid. Journal of Fluid Mechanics. 1965;**22**:709-720. DOI: 10.1017/S002211206500109X

[22] Feng J, Hu HH, Joseph DD. Direct simulation of initial-value problems for the motion of solid bodies in a Newtonian fluid. 1. Sedimentation. Journal of Fluid Mechanics. 1994;**261**: 95-134. DOI: 10.1017/S0022112094 000285

[23] Richardson JF, Zaki WN. Sedimentation and fluidization: Part I. Transactions of the Institution of Chemical Engineers. 1954;**32**:35-53. DOI: 10.1016/S0263-8762(97)80006-8

[24] Guazzelli E, Morris JF. A Physical Introduction to Suspension Dynamics. Cambridge: Cambridge University Press; 2012

## Efficient Simulation Tools (EST) for Sediment-Laden Shallow Flows

*Sergio Martínez-Aranda and Pilar García-Navarro*

#### **Abstract**

Rapid flows of water-sediment mixtures are probably the most challenging and unknown geophysical gravity-driven processes. The fluidized material in motion consists of a mixture of water and multiple solid phases with different specific characteristics. Modeling sediment transport involves an increasing complexity due to the variable bulk properties in the sediment-water mixture, the coupling of physical processes, and the presence of multiple layer phenomena. Twodimensional shallow-type mathematical models are built in the context of free surface flows and are applicable to most of these geophysical surface processes. Their numerical solution in the finite volume framework is governed by the dynamical properties of the equations, the coupling between flow variables and the computational grid. The complexity of the numerical resolution of these highly unsteady flows and the computational cost of simulation tools increase considerably with the refinement of the non-structured spatial discretization, so that the computational effort required is one of the biggest challenges for the application of depth-averaged 2D models to large-scale long-term flows. Throughout this chapter, the combination of 2D mathematical models, robust numerical methods, and efficient computing kernels is addressed to develop Efficient Simulation Tools (EST's) for environmental surface processes involving sediment transport with realistic temporal and spatial scales.

**Keywords:** sediment-laden flows, depth-averaged rheological models, bed-material entrainment, shallow flow models, GPU high-performance-computing

#### **1. Introduction**

Sediment transport is ubiquitous in environmental water bodies such as rivers, floods, coasts and estuaries, but also is the main process in natural ladslides, debris flows, muddy slurries or mining tailings (**Figure 1**). These are considered highly solidladen fluids, where the density of the water-solid mixture can be more than twice or three times the water density and the bulk solid phase represents about 40–80% of the flow volume [1].

The presence of the solid phases, especially the fine material as silt or clay, affects the rheological behavior of the mixture. The water-sediment mixture rheology begins

**Figure 1.** *Geophysical surface flows involving sediment transport.*

to be affected by fine solid particle transported in the flow when the volumetric concentration of fine sediment particles reaches about 4% by volume [2], creating a slight shear strength within the fluid. For higher concentrations, the mixture shows a marked non-Newtonian rheology. Mud and debris flows lie between hyperconcentrated flows and wet avalanches [3]. High concentrations of solids generate a critical yield stress which allows that gravel particles can be suspended indefinitely into the fluidized material.

The mass exchange between mud/debris mixtures and erodible beds involves complicated physical processes and the understanding of its theoretical basis remains unclear. Experiments in large-scale channel [4, 5] and field observations in real debris events [6, 7] indicate that the entrainment volume in steep beds can be in the same order of magnitude as the initial volume mobilized. Debris and mud flows gain much of their mass and momentum as they move over steep slopes as a consequence of the material entrainment from the erodible bed, before deposition begins on flatter terrain downstream.

The mathematical modeling of solid-liquid mixture flows and their numerical resolution is still a challenging topic, especially when dealing with realistic applications. When liquid and solid phases are well-mixed, assuming that the solid phase is distributed uniformly over the flow column allows the use of depthaveraged models derived from the vertical integration of the Navier-Stokes equations [8]. Shallow-type mathematical models represent a simplified formulation, derived from the general 3D Navier-Stokes equations, which is applicable to a large number of these geophysical surface processes involving sediment transport. The simplest models, used in river and coastal dynamics, assume small enough sediment concentrations throughout the flow to consider the bulk density constant and uniform. Most of the numerical models reported for highly solid-laden flows also use this one-single-phase approach, neglecting the bulk density in the shallow-flow mass and momentum equations [9, 10]. Nevertheless, even small density gradients influence importantly the mixing dynamics in flow confluences [11] and larger

gradients can also generate numerical oscillations and instabilities throughout mixing interfaces.

Modeling sediment transport involves an increasing complexity with respect to rigid-bed shallow water models [12] due to the presence of variable sediment-fluid mixture properties, coupling of physical processes and multiple layers phenomena [13]. One of the biggest challenges for the application of depth-averaged models to realistic large-scale long-term flows is the computational effort required. New strategies to reduce the computational effort have been developed in the last decade through the use of parallelization techniques based on Multiprocessing (OpenMP) or Message Passing Interface (MPI), which allow to run simulations on multi-CPU clusters. Their main drawback is the associated hardware cost and energy requirements, which are directly proportional to the number of CPU-cores available and limit their efficiency. In the last years, the usage of Graphics Processing Units (GPU) hardware accelerators for sequential computation has demonstrated to be an efficient and low cost alternative to the traditional multi-CPU strategies [14]. GPU-accelerated algorithms have been developed for real-time floods forecasting [15], real-scale bedload erosive shallow-flows [16] or tsunami prediction [17]. GPU devices are oriented to perform arithmetical operations on vector-based information. Unlike the conventional shared-memory multi-CPU implementations, the GPU solution must be designed taking into account the fact that the GPU is an independent device with its own RAM memory. This means that the memory transfer between the conventional RAM memory and the GPU device memory plays a key role in the performance of GPU-accelerated software.

#### **2. Depth-integrated equations for shallow flows**

#### **2.1 Mass and linear momentum conservation**

The flow of a water-sediment mixture can be mathematically described assuming the movement of the solid particles as a diffusion phenomenon into the liquid phase. Then, the continuity and momentum conservation for the mixture, supplemented with the transport equation for the solid phase, can be established for modeling these two-phase flows. Although both solid and liquid phases are incompressible when considered independently, the bulk behavior of the solid-liquid mixture is the same as that of a compressible material depending on the local solid phase volumetric concentration. It possible to define the bulk density *ρ* ¼ *ρwn* þ *ρsϕ* and linear momentum *ρ***u** ¼ *ρwn***uw** þ *ρsϕ***us** of the mixture, being *ρw* and *ρs* the density of the fluid and solid phases, respectively, *ϕ* the volumetric solid-phase concentration and *n* ¼ 1 � *ϕ* the volumetric fluid-phase fraction or mixture porosity, **uw** the velocity of the pore-fluid and **us** the advective sediment particle velocity.

The 3D time-averaged Navier-Stokes equations for mass and momentum conservation of a two-phase mixture can be written in the Cartesian coordinate system **X** ¼ ð Þ *x*, *y*, *z* as

$$\frac{\partial \rho}{\partial t} + \nabla \cdot (\rho \mathbf{u}) = \mathbf{0} \tag{1}$$

$$\frac{\partial(\rho \mathbf{u})}{\partial t} + \nabla \cdot (\rho \mathbf{u} \otimes \mathbf{u}) = \mathbf{F} - \nabla p + \nabla \cdot \mathbf{\tau} \tag{2}$$

where **u** ¼ *ux*, *uy*, *uz* � � is the bulk velocity in any point of the fluidized material, **F** ¼ *Fx*, *Fy*, *Fz* � � are the external forces, such as gravity and *p* denotes the pressure of the mixture. The term *τ* ¼ *τij* ð Þ *i*, *j* ¼ *x*, *y*, *z* is the deviatoric stress tensor. With low solid-phase concentrations, the water-sediment mixture behaves as a Newtonian fluid, with a constitutive relation given by the Navier-Poisson law [8]. Nevertheless, for high sediment concentrations the mixture becomes a kind of non-Newtonian fluid with a complex constitutive law, which depends on multiple factors, relating stresses and deformation rates.

In order to develop a shallow-type depth-averaged mathematical model, the Navier-Stokes system (1) and (2) is integrated between the free surface *zs* ¼ *zs*ð Þ *t*, *x*, *y* and the bottom surface of the flow column *zb* ¼ *zb*ð Þ *t*, *x*, *y* , which is also considered a movable interface. The kinematic conditions at these boundaries can be expressed as

$$\frac{\partial \mathbf{z}\_s}{\partial t} + (u\_\mathbf{x})\_s \frac{\partial \mathbf{z}\_s}{\partial \mathbf{x}} + \left(u\_\mathbf{y}\right)\_s \frac{\partial \mathbf{z}\_s}{\partial \mathbf{y}} = (u\_\mathbf{z})\_s \tag{3}$$

$$\frac{\partial \mathbf{z}\_b}{\partial t} + (u\_x)\_b \frac{\partial \mathbf{z}\_b}{\partial x} + \left(u\_\mathcal{\mathcal{V}}\right)\_b \frac{\partial \mathbf{z}\_b}{\partial y} = (u\_x)\_b + N\_b \tag{4}$$

being *Nb* the net volumetric flux through the bed interface along the *z*�coordinate, being positive for net entrainment conditions and negative for net deposition fluxes. The subscripts ð Þ� *s* and ð Þ� *b* indicate the value of the corresponding variable at the flow free surface and the bottom bed interface respectively.

The mass conservation Eq. (1) is integrated along the water column as

$$\int\_{x\_b}^{x\_t} \partial\_t(\rho) \, \mathrm{d}z + \int\_{x\_b}^{x\_t} \partial\_x(u\_x \rho) \, \mathrm{d}z + \int\_{x\_b}^{x\_t} \partial\_\gamma(u\_\gamma \rho) \, \mathrm{d}z + \int\_{x\_b}^{x\_t} \partial\_z(u\_x \rho) \, \mathrm{d}z = 0 \tag{5}$$

Applying Leibnitz's rule to each term of (5) and using (3) and (4) leads to

$$\frac{\partial(\overline{\rho}h)}{\partial t} + \frac{\partial}{\partial x}(\overline{\rho}h\overline{u}) + \frac{\partial}{\partial y}(\overline{\rho}h\overline{v}) = -(\rho)\_b N\_b \tag{6}$$

where the depth-averaged of the mixture bulk density *ρ* and velocities are defined as

$$
\overline{\rho}h = \int\_{x\_b}^{z\_t} \rho \,\mathrm{d}z \qquad \overline{u} = \frac{1}{\overline{\rho}h} \int\_{x\_b}^{z\_t} u\_x \rho \,\mathrm{d}z \qquad \overline{v} = \frac{1}{\overline{\rho}h} \int\_{x\_b}^{z\_t} u\_y \rho \,\mathrm{d}z \tag{7}
$$

being *h* ¼ *zs* � *zb* the flow depth and ð Þ*ρ b* the mixture bulk density at the bed interface *zb*.

Regarding the momentum depth integration, the volumetric force components in (1) and (2) are considered for the sake of simplicity as *Fx* ¼ *Fy* ¼ 0 and *Fz* ¼ �*ρg*. Moreover, the z-momentum equation in (2) can be simplified by assuming shallow-flow scaling and neglecting temporal, convective and stress terms, leading to a hydrostatic pressure distribution along the flow column.

The x-momentum in (2) is integrated throughout the flow depth in the following way

*Efficient Simulation Tools (EST) for Sediment-Laden Shallow Flows DOI: http://dx.doi.org/10.5772/intechopen.107345*

$$\begin{aligned} \left\| \int\_{\boldsymbol{x}\_{b}}^{\boldsymbol{z}\_{t}} \partial\_{t} (\rho \boldsymbol{u}\_{\boldsymbol{x}}) \, \mathrm{d}\mathbf{z} + \int\_{\boldsymbol{x}\_{b}}^{\boldsymbol{z}\_{t}} \partial\_{\boldsymbol{x}} (\boldsymbol{u}\_{\boldsymbol{x}} \rho \boldsymbol{u}\_{\boldsymbol{x}}) \, \mathrm{d}\mathbf{z} + \int\_{\boldsymbol{x}\_{b}}^{\boldsymbol{z}\_{t}} \partial\_{\boldsymbol{x}} (\boldsymbol{u}\_{\boldsymbol{y}} \rho \boldsymbol{u}\_{\boldsymbol{x}}) \, \mathrm{d}\mathbf{z} \\ &= - \int\_{\boldsymbol{x}\_{b}}^{\boldsymbol{z}\_{t}} \partial\_{\boldsymbol{x}} (\boldsymbol{p}) \, \mathrm{d}\mathbf{z} + \int\_{\boldsymbol{x}\_{b}}^{\boldsymbol{z}\_{t}} \left( \partial\_{\boldsymbol{x}} \boldsymbol{\tau}\_{\mathrm{xx}} + \partial\_{\boldsymbol{\gamma}} \boldsymbol{\tau}\_{\mathrm{xy}} + \partial\_{\boldsymbol{x}} \boldsymbol{\tau}\_{\mathrm{xx}} \right) \, \mathrm{d}\mathbf{z} \end{aligned} \tag{8}$$

Applying Leibnitz's rule and considering the kinematic boundary conditions at both the free surface (3) and the bed interface (4), the left hand side of (8) is integrated as

$$\frac{\partial(\overline{\rho}h\overline{u})}{\partial t} + \frac{\partial}{\partial \mathbf{x}} \left( \overline{\rho}h\overline{u}^2 - \overline{\rho}hD\_{\mathbf{x}\mathbf{x}} \right) + \frac{\partial}{\partial \mathbf{y}} \left( \overline{\rho}h\overline{u}\overline{w} - \overline{\rho}hD\_{\mathbf{x}\mathbf{y}} \right) + (\rho u\_{\mathbf{x}})\_{\mathbf{b}}N\_{\mathbf{b}} \tag{9}$$

with *Dxx*, *Dxy* � � accounting for the depth-averaged dispersion momentum transport due to the non-uniformity of the vertical velocity profile, and defined as

$$D\_{\infty} = -\frac{1}{\overline{\rho}h} \int\_{x\_b}^{x\_t} \rho (u\_{\infty} - \overline{u})^2 \,\mathrm{d}x \tag{10}$$

$$D\_{\mathbf{x}\mathbf{y}} = -\frac{1}{\overline{\rho}h} \int\_{x\_b}^{x\_t} \rho (u\_{\mathbf{x}} - \overline{u}) (u\_{\mathbf{y}} - \overline{v}) \, \mathrm{d}x \tag{11}$$

On the right hand side of (8), using the hydrostatic pressure distribution and assuming that the vertical density gradient is negligible compared with those along the horizontal plane, as occurs in natural debris and mud slurry flows, the integral of the pressure gradient along the *x*�coordinate can be expressed as

$$-\int\_{x\_b}^{x\_r} \partial\_x(p) \, \mathrm{d}x = -\mathop{\rm g}\limits\_{x\_b}^{x\_r} \left(\frac{\partial}{\partial x} \bigg\|\, \rho \, \mathrm{d}x\right) \, \mathrm{d}x = \underbrace{-\frac{\partial}{\partial x} \left(\frac{1}{2} \mathrm{g} \overline{\rho} h^2\right)}\_{\text{Conversative term}} - \underbrace{-\mathrm{g} \overline{\rho} h \frac{\mathrm{d}x\_b}{\mathrm{d}x}}\_{\text{Bed} - \text{pressure term}}\tag{12}$$

separating the pressure gradient term into a conservative component plus a bedpressure component.

The stress terms on the right hand side of (8) are also integrated along the flow column, leading to the final expresion for the depth-integrated momentum equation along the *x*�coordinate

$$\begin{split} \frac{\partial(\overline{\rho}h\overline{u})}{\partial t} &+ \frac{\partial}{\partial \mathbf{x}} \left( \overline{\rho}h\overline{u}^2 + \frac{1}{2} \mathbf{g} \overline{\rho}h^2 \right) + \frac{\partial}{\partial \mathbf{y}} (\overline{\rho}h\overline{u}\overline{w}) = -\mathbf{g}\overline{\rho}h \frac{\partial \mathbf{z}\_b}{\partial \mathbf{x}} + \mathbf{\tau}\_{\mathbf{x}\mathbf{x}} - \mathbf{\tau}\_{\mathbf{b}\mathbf{x}} \\ &+ \frac{\partial}{\partial \mathbf{x}} \left( \overline{\rho}h(T\_{\mathbf{x}\mathbf{x}} + D\_{\mathbf{x}\mathbf{x}}) \right) + \frac{\partial}{\partial \mathbf{y}} \left( \overline{\rho}h(T\_{\mathbf{x}\mathbf{y}} + D\_{\mathbf{x}\mathbf{y}}) \right) \\ &- (\rho u\_{\mathbf{x}})\_{b}N\_{b} \end{split} \tag{13}$$

and, similarly, for the *y*�coordinate. The boundary stress terms in (13) at the free surface *zs* and the bed interface *zb* are

$$
\pi\_{\rm xx} = \left(\pi\_{\rm xx}\right)\_{\rm s} - \left(\pi\_{\rm xx}\right)\_{\rm s} \frac{\partial \mathbf{z}\_{\rm s}}{\partial \mathbf{x}} - \left(\pi\_{\rm xy}\right)\_{\rm s} \frac{\partial \mathbf{z}\_{\rm s}}{\partial \mathbf{y}} \tag{14}
$$

$$-\tau\_{b\mathbf{x}} = -(\tau\_{\mathbf{x}\mathbf{x}})\_b + (\tau\_{\mathbf{x}\mathbf{x}})\_b \frac{\partial \mathbf{z}\_b}{\partial \mathbf{x}} + \left(\tau\_{\mathbf{x}\mathbf{y}}\right)\_b \frac{\partial \mathbf{z}\_b}{\partial \mathbf{x}} \tag{15}$$

where *τsx* denotes the *x*�coordinate component of the wind action at the free surface, whereas *τbx* represents the *x*�coordinate component of the boundary shear stress at the bed interface, opposing to the flow movement. Furthermore, the depthintegrated stress terms along the flow column *Txx*, *Txy* � � are defined as

$$T\_{\rm xx} = \frac{1}{\overline{\rho}h} \int\_{x\_b}^{x\_t} \mathbf{r}\_{\rm xx} \, \mathbf{d}z \quad T\_{\rm xy} = \frac{1}{\overline{\rho}h} \int\_{x\_b}^{x\_t} \mathbf{r}\_{\rm xy} \, \mathbf{d}z \tag{16}$$

#### **2.2 Depth-integrated solid transport equation**

The solid-phase transport process is governed by

$$\frac{\partial(\rho\_s \phi)}{\partial t} + \nabla \cdot (\rho\_s \phi \mathbf{u}\_\mathbf{s}) = \mathbf{0} \tag{17}$$

where **us** ¼ *usx*, *usy*, *usz* � � is the local velocity of the solid particles. Because the solid-phase velocity is not included in the dependent variables of the fluid dynamic system, and assuming that the sediment particles are incompressible and non-porous, (17) is rewritten as

$$\frac{\partial \phi}{\partial t} + \nabla \cdot (\phi \mathbf{u}) = \nabla \cdot [(\mathbf{u} - \mathbf{u}\_\mathbf{s})\phi] \tag{18}$$

where the term on the right hand side accounts for the drag effects caused by the liquid phase on the advective solid flux. Due to the definition of the bulk linear momentum of the mixture, this term is usually neglected and the approximation **u**≈**us** is actually assumed.

The transport equation (18) must also be integrated throughout the entire flow column defining the total solid phase being transported in the flow column as

$$\overline{\phi}h = \int\_{x\_b} \phi \,\mathrm{d}z \tag{19}$$

where *ϕ* is the dept-averaged volumetric solid concentration. Therefore, applying the Leibnitz's integration rule to (18) and considering the kinematic boundary conditions at both the free surface (3) and the bed interface (4), the depth-integrated equation for the solid phase in the flow column reduces to

$$\frac{\partial(\overline{\phi}h)}{\partial t} + \frac{\partial}{\partial \mathbf{x}}(\overline{\phi}h\overline{u}) + \frac{\partial}{\partial \mathbf{y}}(\overline{\phi}h\overline{v}) = -(\phi)\_b N\_b + \frac{\partial}{\partial \mathbf{x}}(\overline{\phi}hD\_{\mathbf{x}}) + \frac{\partial}{\partial \mathbf{y}}(\overline{\phi}hD\_{\mathbf{y}}) \tag{20}$$

where ð Þ *ϕ b* denotes the solid concentration in the bed layer surface *zb* and *Dsx*, *Dsy* � � account for the depth-averaged dispersive flux due to the non-uniformity of both the velocity and solid concentration profiles throughout the flow column, defined as

$$D\_{\infty} = -\frac{1}{\overline{\phi}h} \int\_{x\_b}^{x\_t} \rho (u\_{\infty} - \overline{u}) \left( \frac{\phi}{\rho} - \frac{\overline{\phi}}{\overline{\rho}} \right) \, \mathrm{d}x \tag{21}$$

$$D\_{\mathcal{Y}} = -\frac{1}{\overline{\phi}h} \int\_{x\_b}^{x\_t} \rho \left( u\_{\mathcal{Y}} - \overline{v} \right) \left( \frac{\phi}{\rho} - \frac{\overline{\phi}}{\overline{\rho}} \right) \, \mathrm{d}x \tag{22}$$

Furthermore, assuming that the volumetric solid concentration in the bed layer is 1 � *ξ* with *ξ* denoting the bulk porosity of the bed layer, the solid mass conservation in the bed layer requires that

$$(\mathbf{1} - \boldsymbol{\xi})\frac{\partial \mathbf{z}\_b}{\partial t} = (\boldsymbol{\phi})\_b \mathbf{N}\_b \tag{23}$$

#### **2.3 Constitutive model for complex viscoplastic flows**

So far, there is not a universal closure relation for representing the viscous terms in complex non-Newtonian flows. Constitutive formulations used for environmental sediment-water mixtures are mainly derived from 3D general rheological models which, assuming isotropic material and isochoric flow, allow to express the deviatoric component of the stress tensor *σ* in the material as

$$
\pi = \Phi\_1(I\_{2D}) \mathbf{D} \tag{24}
$$

where **D** � *Dij* ¼ 1 2 *∂jui* þ *∂iuj* � � ð Þ *i*, *j* ¼ *x*, *y*, *z* is the rate of deformation tensor and Φ1 is a scalar function of the second invariant *I*2*D* ¼ 1 2 tr **D**2 � � of the rate of deformation tensor **D** [18] and depends on multiple factor, such as cohesive stress, pore-fluid pressure or flow initial regime. The generalized viscoplastic model, also called Herschel-Bulkley model, assumes dependence of Φ1 on the three parameters: *τ*0 the cohesive-frictional yield strength, *K* a consistency coefficient, and *m* a parameter characterizing the rheological response of the mixture [13, 18], known as behavior index. It is worth mentioning that the dimensions of consistency coefficient *K* depends on the behavior index *m*. Therefore, the function Φ1 is expressed as

$$\Phi\_1(I\_{2D}) = \frac{\tau\_0}{\sqrt{I\_{2D}}} + 2K \left(4I\_{2D}\right)^{\frac{m-1}{2}} \tag{25}$$

Considering simple shear state along the flow direction, the velocity vector throughout the flow column is expressed as **u** ¼ *U z*ð Þ*nux*, *U z*ð Þ*nuy*, 0 � �, where *U z*ð Þ is the modulus of the bulk mixture velocity and **nu** ¼ *nux*, *nuy* � � is the velocity unit vector. Therefore

$$\mathbf{D} = \begin{pmatrix} 0 & 0 & \frac{1}{2}\frac{\mathrm{d}U}{\mathrm{d}z}n\_{\mathrm{ux}} \\ 0 & 0 & \frac{1}{2}\frac{\mathrm{d}U}{\mathrm{d}z}n\_{\mathrm{uy}} \\ \frac{1}{2}\frac{\mathrm{d}U}{\mathrm{d}z}n\_{\mathrm{ux}} & \frac{1}{2}\frac{\mathrm{d}U}{\mathrm{d}z}n\_{\mathrm{uy}} & 0 \end{pmatrix} \quad I\_{2D} = \frac{1}{4}\left(\frac{\mathrm{d}U}{\mathrm{d}z}\right)^{2} \tag{26}$$

Replacing (26) and (25) into (24) leads to

$$
\tau = \begin{pmatrix}
\mathbf{0} & \mathbf{0} & \tau(\mathbf{z})\boldsymbol{n}\_{\text{ux}} \\
\mathbf{0} & \mathbf{0} & \tau(\mathbf{z})\boldsymbol{n}\_{\text{uy}} \\
\tau(\mathbf{z})\boldsymbol{n}\_{\text{ux}} & \tau(\mathbf{z})\boldsymbol{n}\_{\text{uy}} & \mathbf{0}
\end{pmatrix} \tag{27}
$$

being *τ*ð Þ*z* the shear stress along the flow direction, which depends on the fluid rheology function (25). For the Newtonian constitutive model, the yield strength *τ*0 is null, *K Pa* ½ � � *s* is the dynamic viscosity of the fluid *μ* and *m* ¼ 1. The generalized model reduces to *τ*ð Þ¼ *z μ* d*U* d*z* for viscous Newtonian fluids. A widespread non-Newtonian constitutive relation for geophysical surface flows in laminar regime is the Bingham model, which considers a pure cohesive yield stress *τ*0 ¼ *τy* for the flow initiation, *K* ¼ *μB* ½ � *Pa* � *s* the Bingham plastic viscosity and *m* ¼ 1. In this case, the generalized model reduces to *τ*ð Þ¼ *z τy* þ *μB* d*U* d*z* .

The frictional non-linear viscoplastic model considers a Coulomb-Terzaghi linear relation between the effective normal stress *σe*ð Þ*z* and the shear stress, hence

$$\pi\_0 = \sigma\_\epsilon(\mathbf{z}) \tan \delta\_\! \!/ = \left[ \rho \mathbf{g}(\mathbf{z}\_\sharp - \mathbf{z}) - \mathcal{P}(\mathbf{z}) \right] \tan \delta\_\! \!/ \tag{28}$$

where Pð Þ*z* denotes de pore-fluid pressure and *δf* accounts for the effective friction angle bewteen solid particles [19]. Estimation of the pore pressure distribution Pð Þ*z* is a challenging task [1, 19, 20], although its effects on the reduction of the intergranular shear stress seem to be demonstrated [5, 21, 22]. The simplest models divide the pore pressure into a hydrostatic component plus a dynamic additive component Pð Þ¼ *z* ð Þ 1 þ E *ρwg z*ð Þ *s* � *z* , being *ρw* the density of the pore-liquid and E a tunning coefficient which usually takes values from about 0.4 to 0.8.

Using (28), and considering a plastic viscosity *K* ¼ *μP Pa* � *s m* ½ �, the generalized model reduces to

$$\tau(\mathbf{z}) = \sigma\_{\mathbf{c}}(\mathbf{z}) \tan \delta\_f + \mu\_P \left(\frac{\mathbf{d}U}{\mathbf{d}z}\right)^m \tag{29}$$

Formulation of closure models for the bottom shear stress requires to integrate the deviatoric stress tensor (24) throughout the flow column. This is not a trivial problem since the structure of the flow along the vertical direction is lost and only the averaged quantities are available. Assuming the deviatoric stress tensor (27), all the depthaveraged stress terms *Txx*, *Txy* � � and *Tyx*, *Tyy* � � (16) are null, and

$$
\tau\_{b\mathbf{x}} = (\tau\_{\mathbf{x}x})\_b - (\tau\_{\mathbf{x}x})\_b \frac{\partial \mathbf{z}\_b}{\partial \mathbf{x}} - \left(\tau\_{\mathbf{x}y}\right)\_b \frac{\partial \mathbf{z}\_b}{\partial y} = \tau\_b n\_{\mathbf{u}\mathbf{x}} \tag{30}
$$

$$
\tau\_{by} = \left(\tau\_{\rm yx}\right)\_b - \left(\tau\_{\rm yx}\right)\_b \frac{\partial \mathbf{z}\_b}{\partial \mathbf{x}} - \left(\tau\_{\rm yy}\right)\_b \frac{\partial \mathbf{z}\_b}{\partial \mathbf{y}} = \tau\_b \cdot \mathbf{n}\_{u\mathbf{y}} \tag{31}
$$

where *τb* ¼ *τ*ð Þ *zb* is the shear stress at the basal interface along the flow direction. In order to obtain a depth-averaged formulation for the basal shear stress *τb*, the distributed shear stress function *τ*ð Þ*z* must be integrated along the flow column. This allows to obtain both the velocity distribution along the vertical direction and the basal resistance *τb* expressed as a function of the depth-averaged flow variables. If the pore-pressure excess in (28) is considered linear, the the constitutive eq. (29) can be rewritten as

$$\pi(\mathbf{z}) = \tau\_f \left( \mathbf{1} - \frac{\mathbf{z} - \mathbf{z}\_b}{h} \right) + \mu\_P \left( \frac{\mathbf{d}U}{\mathbf{d}z} \right)^m \tag{32}$$

being *τf* ¼ ð Þ *ρgh* � P*b* tan *δf* the value of the frictional yield stress at the basal surface and P*b* is the pore-pressure at the bed surface (**Figure 2**).

Therefore, assuming the induced shear distribution along the flow column follows

$$
\pi(z) = \tau\_b \left( 1 - \frac{z - z\_b}{h} \right) \tag{33}
$$

the velocity derivative along the vertical direction can expressed as

$$\frac{\mathbf{d}U}{\mathbf{d}z} = \left[\frac{\tau\_b - \tau\_f}{\mu\_P} \left(\mathbf{1} - \frac{z - z\_b}{h}\right)\right]^{1/m} \tag{34}$$

and integrating (34) leads to the velocity vertical distribution

$$U(z) = \frac{m}{m+1} \left(\frac{\tau\_b - \tau\_f}{\mu\_P}\right)^{1/m} h \left[1 - \left(1 - \frac{z - z\_b}{h}\right)^{\frac{m+1}{m}}\right] \tag{35}$$

for the non-linear viscoplastic model. Note that the velocity at the free surface *U z*ð Þ� *s Uh* and the depth averaged velocity *U* can be expressed as

$$U\_h = \frac{m}{m+1} \left(\frac{\tau\_b - \tau\_f}{\mu\_P}\right)^{1/m} h \qquad \overline{U} = \frac{m+1}{2m+1} U\_h \tag{36}$$

Integrating (35) throughout the flow column allows to obtain the basal shear stress *τb* as a function of the flow depth *h*, the averaged flow velocity *U* ¼ ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi *u*2 þ *v*2 p and the basal frictional yield strength *τf* using

$$
\pi\_b = \pi\_f + \left(\frac{2m+1}{m}\right)^m \mu\_P \left(\frac{\sqrt{\overline{u}^2 + \overline{v}^2}}{h}\right)^m \tag{37}
$$

**Figure 2.** *Velocity and stress distribution for the generalized non-linear frictional model.*

**Figure 3.** *Basal resistance behavior for the generalized frictional non-linear viscoplastic model.*

It is worth mentioning that (37) represents a generalized depth-integrated formulation for viscoplastic flows (**Figure 3**) which encompasses: pseudoplastic or shearthinning behavior for *m* <1, reducing the apparent viscosity as the induced shear rate increases; linear viscoplastic behavior for *m* ¼ 1, with a linear relation between shear stress and shear rate; and dilatant or shear-thickening behavior for *m* >1, increasing the apparent viscosity as the induced shear rate grows.

#### **2.4 Net mass exchange between bed and flow layers**

The net volumetric exchange *Nb* between the underlying bed layer and the mixture flow column at the bed surface *zb*ð Þ *t*, *x*, *y* , appearing in (6), (20) and (23), is usually modeled as the balance between the entrainment and the deposition vertical solid fluxes, *Eb* and *Db* respectively, leading to

$$N\_b = \frac{1}{(\phi)\_b} (D\_b - E\_b) \tag{38}$$

involving the bed surface solid concentration ð Þ *ϕ b* ¼ 1 � *ξ*.

The deposition rate *Db* is commonly related to the solid particles settling velocity in the mixture *ωsm* and to the near-bed solid concentration in the flow *ϕz*!*zb* , which is usually related to the depth-averaged solid concentration in the flow as *ϕz*!*zb* ¼ *αϕ*, being *α* an adaptation or recovery coefficient, leading to *Db* ¼ *αωsm ϕ*.

The settling velocity of the solid particles *ωsm* in highly concentrated mixtures is influenced by the presence of other solid particles. Furthermore, in dense-packed mixtures with moderate plastic fine fractions in the flow column such as muddy slurries, the particle settling velocity can be strongly reduced by the development of internal yield stresses in the pore-fluid. Richardson and Zaki [23] proposed *ωsm* ¼ 1 � *ϕ mωs*, being *ωs* the settling velocity of the sediment particles in clear water and *m* a hindering empirical exponent depending on the Reynolds particle number ( *Re p* ¼ *ωs ds=ν*, with *ν* the clear water kinematic viscosity) which usually takes values close to *m* ¼ 4. Therefore *Efficient Simulation Tools (EST) for Sediment-Laden Shallow Flows DOI: http://dx.doi.org/10.5772/intechopen.107345*

$$D\_b = a\_\mathcal{D} a\_\mathfrak{s} \overline{\phi} \left(\mathbf{1} - \overline{\phi}\right)^4 \tag{39}$$

where *αD* is a dimensionless parameter which accounts for mulple factor and requires calibration.

The erosion solid flux *Eb* is directly related to the turbulent fluctuation of the volumetric solid concentration and flow velocity near the bed surface. We assume that this near-bed erosion rate is at the capacity of the flow to entrain solid material from the underlying bed layer, hence it is related to the settling velocity of the particles in clear water *ωs* and the near equilibrium concentration *ϕ*∗ *z*!*zb* . The near-bed equilibrium concentration is related to the depth-averaged equilibrium concentration *ϕ*∗ as *ϕ*∗ *z*!*zb* ¼ *α*∗ *ϕ*∗ , being *α*∗ an adaptation coefficient under equilibrium conditions. When equilibrium solid transport states are reached, the adaptation coefficients *α* and *α*∗ coincide but in non-equilibrium states *α*∗ 6¼ *α* generally. Therefore, for the sake of simplicity, we assume that the vertical erosion rate *Eb* is expressed as

$$E\_b = \alpha\_{\rm E} \alpha\_{\rm s} \overline{\phi^\*} \tag{40}$$

where *αE* is a dimensionless empirical parameter which requires calibration. The capacity solid concentration is usually computed as *ϕ*∗ ¼ *q* ∗ *s = h* ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi *u*2 þ *v*2 � � p , where *q* ∗ *s* accounts for the value of the solid transport throughout the flow column in capacity or equilibrium condition, which can be estimated using the multiple empirical relationships from the local hydrodynamic variables [8].

#### **3. Efficient numerical tools for multi-grain mud/debris flows**

Considering a multi-grain mixture flow, the resulting system is composed by 3 þ *N* þ 1 conservation equations accounting for the bulk mass (6) and momentum (13), the transport of the *N* sediment classes (20) and the bed elevation evolution (23). The dimensionless bulk density *r* can be expressed by defining a new variable *ϕχ* , referred to as buoyant solid concentration

$$r = \frac{\rho}{\rho\_w} = \mathbf{1} + \phi^\chi \quad \text{with} \quad \phi^\chi = \sum\_{p=1}^N \frac{\rho\_{s,p} - \rho\_w}{\rho\_w} \phi\_p \tag{41}$$

where *ρs*,*p* and *ϕp* are the density and depth-averaged volumetric concentration of the *p*th solid phase respectively. Using (41), the equations forming the system can be recast as five conservation laws and rewritten in vector form as

$$\frac{\partial \mathbf{U}}{\partial t} + \nabla \cdot \mathbf{E}(\mathbf{U}) = \mathbf{S}\_{\mathbf{b}}(\mathbf{U}) + \mathbf{S}\_{\mathbf{r}}(\mathbf{U}) + \mathbf{E}\_{\mathbf{b}}(\mathbf{U}) \tag{42}$$

where **U** is the vector of conservative variables and **E U**ð Þ¼ ð Þ **F U**ð Þ, **G U**ð Þ are the convective fluxes along the **X** ¼ ð Þ *x*, *y* horizontal coordinates, respectively, expressed as

$$\mathbf{U} = \begin{pmatrix} rh \\ rhu \\ rhv \\ h\phi^{\mathbb{X}} \\ z\_b \end{pmatrix} \quad \mathbf{F}(\mathbf{U}) = \begin{pmatrix} rhu \\ rhu^2 + \frac{1}{2}grh^2 \\ rhuv \\ rhu\phi^{\mathbb{X}} \\ 0 \end{pmatrix} \quad \mathbf{G}(\mathbf{U}) = \begin{pmatrix} rhv \\ rhuv \\ rhv^2 + \frac{1}{2}grh^2 \\ hv\phi^{\mathbb{X}} \\ 0 \end{pmatrix} \tag{43}$$

The vector **Sb**ð Þ **U** accounts for the momentum source term associated to the variation of the pressure force on the bed interface, whereas **S***τ*ð Þ **U** is the momentum dissipation due to the basal resistance. Finally, the source term **Eb**ð Þ **U** accounts for the bulk mass exchange between the mixture flow and the bed layer.

$$\mathbf{S}\_{\mathbf{b}}(\mathbf{U}) = \begin{pmatrix} 0\\ -g\eta h \frac{d\mathbf{z}\_{b}}{\alpha\varepsilon} \\ -g\eta h \frac{d\mathbf{z}\_{b}}{\partial\mathbf{y}} \\ 0 \\ 0 \\ 0 \end{pmatrix}, \quad \mathbf{S}\_{\mathbf{r}}(\mathbf{U}) = \begin{pmatrix} 0\\ -\frac{\tau\_{b}}{\rho\_{w}} \eta\_{w\mathbf{x}} \\\ -\frac{\tau\_{b}}{\rho\_{w}} \eta\_{w\mathbf{y}} \\\ 0 \\\ 0 \\ 0 \end{pmatrix}, \quad \mathbf{E}\_{\mathbf{b}}(\mathbf{U}) = \begin{pmatrix} -\frac{\rho\_{b}}{\rho\_{w}(1-\xi)} \sum\_{p=1}^{N} \left(D\_{b} - E\_{b}\right)\_{p} \\\ 0 \\\ 0 \\\ -\sum\_{p=1}^{N} \frac{\rho\_{s,p} - \rho\_{w}}{\rho\_{w}} \left(D\_{b} - E\_{b}\right)\_{p} \\\ -\frac{1}{\tau} \sum\_{p=1}^{N} \left(D\_{b} - E\_{b}\right)\_{p} \\\ \frac{1}{\mathbf{1}-\xi} \sum\_{p=1}^{N} \left(D\_{b} - E\_{b}\right)\_{p} \end{pmatrix} \tag{44}$$

#### **3.1 Finite volume method in unstructured meshes**

System (42) is time dependent, non linear and contains mass and momentum source terms. Under the hypothesis of dominant advection it can be classified as belonging to the family of hyperbolic systems. In order to obtain a numerical solution, the spatial domain is divided in computational cells using a fixed-in-time mesh and system (42) is integrated in each cell Ω*i*. Applying the Gauss theorem leads to

$$\frac{\partial}{\partial t} \int\_{\Omega\_i} \mathbf{U} \, \mathbf{d} \Omega + \sum\_{k=1}^{\text{NE}} (\mathbf{E} \cdot \mathbf{n})\_k \, l\_k = \int\_{\Omega\_i} \mathbf{S}\_\mathbf{b}(\mathbf{U}) \, \text{d}\Omega + \int\_{\Omega\_i} \mathbf{S}\_\tau(\mathbf{U}) \, \text{d}\Omega + \int\_{\Omega\_i} \mathbf{E}\_\mathbf{b}(\mathbf{U}) \, \text{d}\Omega \tag{45}$$

being *NE* the number of edges for the *i* cell, **n** ¼ *nx*, *ny* � � *k* the outward unit normal vector, *lk* the length of the edge and hence ð Þ **E** � **n** *k* the value of the normal flux through the *k*th edge (**Figure 4**).

The left hand side of (42), also called homogeneous part, satisfies the rotation invariant property [24] and hence the conservative normal flux vector can be rewritten as

$$\mathbf{(E \cdot n)}\_{k} = \left[ \mathbf{F(U)} n\_{\rm x} + \mathbf{G(U)} n\_{\rm y} \right]\_{k} = \mathbf{R}\_{k}^{-1} \mathbf{F(R\_{k} U)} \tag{46}$$

being **R** a rotation matrix which projects the global orthogonal framework **X** ¼ ð Þ *x*, *y* into the local framework **X**^ ¼ ð Þ *x*^, ^*y* of the *k*th cell edge (**Figure 4**),

*Efficient Simulation Tools (EST) for Sediment-Laden Shallow Flows DOI: http://dx.doi.org/10.5772/intechopen.107345*

**Figure 4.** *Computational cells in triangular meshes and local coordinates at the kth cell edge.*

corresponding to the normal and the tangential directions to the edge respectively. The rotation matrix **R***k* and its inverse **R**�**1** *k* are defined as

$$\mathbf{R}\_{k} = \begin{pmatrix} 1 & 0 & 0 & 0 & 0 \\ 0 & n\_{\times} & n\_{\mathcal{V}} & 0 & 0 \\ 0 & -n\_{\mathcal{V}} & n\_{\times} & 0 & 0 \\ 0 & 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 0 & 1 \end{pmatrix}\_{k} \mathbf{R}\_{k}^{-1} = \begin{pmatrix} 1 & 0 & 0 & 0 & 0 \\ 0 & n\_{\times} & -n\_{\mathcal{V}} & 0 & 0 \\ 0 & n\_{\times} & n\_{\times} & 0 & 0 \\ 0 & 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 0 & 1 \end{pmatrix}\_{k} \tag{47}$$

and the set of local conservative variables **U**^ and fluxes **F U**^� � is defined as

$$\hat{\mathbf{U}} \equiv \mathbf{R}\_k \mathbf{U} = \begin{pmatrix} rh \\ rh\hat{u} \\ rh\hat{v} \\ h\phi^\mathbf{r} \\ z\_b \end{pmatrix} \qquad \mathbf{F}(\hat{\mathbf{U}})\_k \equiv \mathbf{F}(\mathbf{R}\_k \mathbf{U}) = \begin{pmatrix} rh\hat{u} \\ rh\hat{u}^2 + \frac{1}{2}\mathbf{g}\_n rh^2 \\ rh\hat{u}\hat{v} \\ h\hat{u}\phi^\mathbf{r} \\ \mathbf{0} \end{pmatrix} \tag{48}$$

where *u*^ ¼ *unx* þ *vny* and ^*v* ¼ �*uny* þ *vnx* are the components of the flow velocity **u**^ in the local framework.

The value of the fluxes through the *k*th cell edge can be augmented incorporating the contribution of the momentum source terms **Sb** and **S***τ* into the homogeneous normal fluxes **F U**^� � [25]. The momentum terms can be included within the local framework ð Þ *x*^, ^*y* using the spatial discretization

$$\int\_{\Omega} \mathbf{S}\_{\mathbf{b}}(\mathbf{U}) \, \mathrm{d}\Omega = \sum\_{k=1}^{\mathrm{NE}} \mathbf{R}\_{k}^{-1} \mathbf{H}(\hat{\mathbf{U}})\_{k} \, l\_{k} \qquad \int\_{\Omega} \mathbf{S}\_{\mathbf{r}}(\mathbf{U}) \, \mathrm{d}\Omega = \sum\_{k=1}^{\mathrm{NE}} \mathbf{R}\_{k}^{-1} \mathbf{T}(\hat{\mathbf{U}})\_{k} \, l\_{k} \tag{49}$$

where **H U**^� � *k* and **T U**^� � *k* are the integrated bed pressure [26] and basal resistance [27] throughout the *k*th cell edge, expressed in the local framework, being

$$\mathbf{H}(\hat{\mathbf{U}})\_k = \begin{pmatrix} \mathbf{0}, \ -\operatorname{gph}\Delta\mathbf{z}\_b, \ \mathbf{0}, \ \mathbf{0}, \ \mathbf{0} \end{pmatrix}\_k^T \tag{50}$$

$$\mathbf{T(\hat{U})}\_k = \left(\mathbf{0}, -\frac{\tau\_b}{\rho\_w} \left(n\_{w\mathbf{x}} \Delta \mathbf{x} + n\_{w\mathbf{y}} \Delta \mathbf{y}\right), \mathbf{0}, \mathbf{0}, \mathbf{0}\right)\_k^T \tag{51}$$

Using (46) and (49), the augmented flux at the *k*th edge can be defined as

$$\mathcal{F}(\hat{\mathbf{U}})\_k = \left(\mathbf{F}(\hat{\mathbf{U}}) - \mathbf{H}(\hat{\mathbf{U}}) - \mathbf{T}(\hat{\mathbf{U}})\right)\_k \tag{52}$$

The net exchange flux term **Eb**ð Þ **U** is discretized in space as

$$\int\_{\Omega\_i} \mathbf{E}\_{\mathbf{b}}(\mathbf{U}) d\Omega \approx A\_i \mathbf{E}\_{\mathbf{b}}(\mathbf{U}\_i) \tag{53}$$

and, assuming a piecewise uniform representation of the variables at the *i* cell for the time *t* ¼ *t n* and using explicit temporal integration for the mass and momentum source terms, the updating formulation for the conservative variables **U** is expressed as

$$\mathbf{U}\_{i}^{n+1} = \mathbf{U}\_{i}^{n} - \frac{\Delta t}{A\_{i}} \sum\_{k=1}^{\text{NE}} \mathbf{R}\_{k}^{-1} \mathcal{F} \left( \hat{\mathbf{U}}\_{i}^{n}, \hat{\mathbf{U}}\_{j}^{n} \right)\_{k}^{\downarrow} l\_{k} + \Delta t \mathbf{E}\_{\mathbf{b}} \left( \mathbf{U}\_{i}^{n} \right) \tag{54}$$

being Δ*t* ¼ *t n*þ1 � *t n* the time step, *Ai* the discrete cell area and F **U**^ *n i* , **U**^ *n j* � �↓ *k* the

numerical flux through the *k*th cell edge computed as a function of the local conservative variables at the neighboring cells *i* and *j*. The numerical fluxes are here upwind computed using a fully-coupled Roe-type Riemann solver (RS) adapted to compressible two-phase shallow flows. Note that the flow density and depth remain coupled in both the conservative variables and fluxes on the left hand side of the equations, improving the robustness and accuracy of the solution when large density gradients appear in the flow [28]. The RS formulation is based on the augmented Roe strategy [25], ensuring the well-balance in steady states and the correct treatment of wet-dry fronts without requiring additional time step restrictions. A detailled explanation of this fully-coupled RS can be found in Aranda et al. [29].

#### **3.2 High-performance-computing in GPU's**

Graphics Processing Units (GPU's), were developed to control and manage the graphic operations for video games but currently they have become a powerful tool for solving engineering problems. The main advantage resides in the high number of computation nodes available for workload distribution, leading to higher speed-ups without an increment of investment in large facilities. The main drawback is that GPU-kernels usually require an intensive programming effort compared with multi-CPU strategies, such as OpenMP or other shared-memory strategies. There exist several platforms for High-Performance-Computing (HPC). NVIDIA developed the CUDA (Compute Unified Device Architecture), a computing platform and programming environment which incorporates directives for massive parallel operations. To ensure the applicability of the model to realistic large-scale mud/debris flows, the updating eq. (54) is solved using a GPU-based algorithm implemented in the NVIDIA CUDA/C++ framework.

The preprocess step and CPU-GPU memory transfer are implemented to run on one CPU core, whereas the time loop computation is accelerated using GPU. However, some tasks inside the time loop are controlled yet by the CPU, such as the time advance control, the boundary conditions application and the output data dump. Therefore, it is necessary to transfer information from/to the GPU at each time-step. *Efficient Simulation Tools (EST) for Sediment-Laden Shallow Flows DOI: http://dx.doi.org/10.5772/intechopen.107345*

While the computational effort required for the time and boundaries transference is considerably smaller than that of each kernel function, in order to dump the intermediate output information, all the variables in the domain must be transferred from GPU device to CPU host.

The CUDA toolkit allows that all the processed elements can be distributed by threads and blocks of threads. Each thread uses its own thread index to identify the element to be processed, launching several execution threads at the same time (parallel computation). As computing GPU devices are well designed to work efficiently with ordered information, the variables needed for computation are stored in the GPU memory as structures of arrays (SoA), improving the spatial locality for memory accesses. Only the kernel functions, which require a higher computational effort, have been implemented to run on the GPU device. Some tasks in the GPU kernel are optimized using the CUBLAS library included in CUDA. The memory transfer between the CPU host and the GPU device has been reduced as much as possible for each time step.

#### **4. Application to the mine tailings dam failure in Brumadinho (Brasil)**

In this section, the proposed methodology is applied to the catastrophic large-scale mud flow occurred Brumadinho (Minas Gerais, Brazil, 2019). The sudden failure of a mine containing dike resul in an extremely violent tailings flow which traveled downstream more than 10 km and reached the Paraopeba River. This disaster caused more than 260 deaths and important economic and environmental losses. The dam contained almost 12 106 *m*3 mining waste tailings with a height of 70 80*m* and covering an area of 4*:*13 105 *m*2. **Figure 5** shows an aerial image of the mine site after

**Figure 5.** *Aerial image of the area affected by the mud.*

#### **Table 1.**

*Brumadinho's dam and tailings features.*

the dam collapse. The thalweg elevation along the area covered by the mud varies between 860*m:o:s:l* at the dam-toe and 720 *m:o:s:l* at the Paraopeba River, with an averaged longitudinal bed slope *S*0 ¼ 0*:*0165 *m=m*. The area affected by the mud was 3*:*3 � 106 *m*2.

The dike material and the tailings rapidly became a heavy liquid that flowed downstream at a high speed. Tailings were composed by a mixture of water, sediments and heavy metals, mainly iron, aluminum, manganese and titanium [30]. The size distribution consisted basically of a mineral sand fraction (38%) and a fines fraction (62%), accounting for mineral silt-clay and metals particles. The water content before the failure was estimated around 50% by volume with a specific weight of 22 � 26 kN*=*m3 (**Table 1**).

The spatial domain is discretized using a unstructured triangular mesh with 6 � 105 cells approximately and four control cross-sections are placed downstream of the dam (**Figure 5**), corresponding to (CS-1) the mine stockpile area, (CS-2) the railway bridge, (CS-3) the national road and (CS-4) the gauge station in the river. Furthermore, the base regime water depth and velocities in the Paraopeba River before to the mud arrival is assessed by a previous simulation. The initial tailings depth at the dam is estimated by comparing the terrain elevation before and after the dam failure using 1 � 1 m DTM's. Six different solid phases are set for characterizing the fully saturated tailings material, including mineral sand, mineral silt, iron (Fe), aluminum (Al), manganese (Mn) and titanium (Ti), with initial bulk volumetric concentration *ϕ*0 ¼ 0*:*5 and normalized density *r*≈2*:*25. The deposition porosity *ξ* for the bed layer is estimated using the Wu relation [8] and the hiding-exposure effects on the critical Shields stress for the incipient motion of each solid phase *θc*,*p* are estimated using the Egiazaroff formula [31]. Also, the solid transport capacity of the flow *q* ∗ *s* is computed using the Wu relation [8]. The simulated time was 3 h from the dam collapse.

The mining tailings in the dam showed a low plasticity and high values of porefluid pressure [30] before the collapse. Assuming the turbulent behavior of the flow, the behavior index is set *m* ¼ 2, setting a quadratic relation in (37) and allowing to calculate the apparent consistency parameter *μP* as a function of the Manning's roughness coefficient of the terrain (*nb* ¼ 0*:*065*sm*�1*=*3). The basal stability angle is estimated in *δf* ¼ 5° . The basal pore pressure is estimated as P*b* ¼ ð Þ 1 þ E*b ρwgh*, being E*b* a basal pressure factor which tends to E*b* ! ∞ for totally fluidized materials. Therefore, the total frictional-turbulent basal shear stress is estimated as

*Efficient Simulation Tools (EST) for Sediment-Laden Shallow Flows DOI: http://dx.doi.org/10.5772/intechopen.107345*

$$\frac{\tau\_b}{\rho\_w} = (r - 1 - \mathcal{E}\_b)gh + rg \frac{n\_b^2}{h^{1/3}} \left( u^2 + v^2 \right) \tag{55}$$

where the first term on the right hand side accounts for the frictional yield stress whereas the second term denotes the turbulent shear component.

In order to assess the influence of the basal pressure in the numerical results, the entrainment term **Eb**ð Þ **U** in (42) is neglected and the value of E*b* is varied from E*b* ! ∞ (pure Newtonian turbulent behavior) to E*b* ¼ 0*:*9 (medium frictional yield strength). **Figure 6** shows the temporal evolution of the wave-front location and the total flow rate at the control section CS-2. Note that, for each value of E*b*, the wave front loses velocity progressively as it moves downstream. As E*b* decreases, the frictional shear stress is enhanced and the mobility of the flow decreases considerably. That it is clearly shown in the marked reduction of the flow rate at the control section CS-2 as the pore pressure factor is increased. Setting E*b* ¼ 1*:*05 allows to predict the observed arrival time (32 � 47 min ) of the tailings wave to the Paraopeba river at 8*:*5 km downstream the dam (gray rectangle in **Figure 6**).

**Figure 7** shows the mud depth at *t* ¼ 0 min , *t* ¼ 5 min , *t* ¼ 15 min and *t* ¼ 45 min after the dam collapse for E*b* ¼ 1*:*05. The numerical results show that practically the whole initial tailing volume flows out of the dam in the first 5 min after the dam collapse, as it was observed in the available videos. The mud wave moves downstream with a computed height larger than 20 m in some zones during the first minutes. After this initial stage, the numerical results show that the dambreak wave is still higher than 10 m when the wave-front reaches the railway bridge (CS-2). During the real event, the mud wave impact caused the collapse of the railway bridge structure. At *t* ¼ 45 min , the mud wave-front has reached the Paraopeba River and the flow is practically stopped.

Furthermore, in order to assess the influence of the solid material entrainment from the erodible bed to the flow layer, the pore-pressure factor is set to E*b* ¼ 1*:*05 but the entrainment term **Eb**ð Þ **U** in (42) is now considered. The relation *αE=αD*, which controls the balance between the erosion component *Eb* and the deposition component *Db* of the vertical bed-flow exchange flux, is varied to *αE=αD* ¼ 1%, *αE=αD* ¼ 5%, *αE=αD* ¼ 10% and *αE=αD* ¼ 50%. **Figure 8** depicts the predicted temporal evolution of both the mud wave-front location and the total volume of the fluidized material in the domain as the entrainment parameter is varied. The higher the relation *αE=αD*, the larger the fluidized mass involved in the flow as the dambreak wave progresses. Hence, the entrainment of bed material leads to a greater flow mobility, causing the wave-front to reach the Paraopeba River faster as the entrainment parameter *αE=αD*

**Figure 6.** *Temporal evolution of (left) the wave-front location and (right) the total flow rate at the control section CS-2.*

**Figure 7.** *Mud depth at t* ¼ 0 min *, t* ¼ 5 min *, t* ¼ 15 min *and t* ¼ 45 min *after the dam collapse.*

#### **Figure 8.**

*Temporal evolution of (left) the wave-front location and (right) the total fluidized volume as the relation αE=αD is increased.*

increases. Even, for *αE=αD* ¼ 50%, the mobilized mass becomes more than 150% of the initial mass and the wave-front reaches the river faster than the fixed-bed pureturbulent case (dashed gray line).

Furthermore, the entrainment of bed materials also leads to the segregation of the solid phases along the flow. As *αE=αD* is increased, horizontal gradients in the concentration of the different solid classes composing the mixture tend to appear (**Figure 9**). These gradients depends on the temporal evolution of the class-specific erosion and deposition vertical fluxes, *Ep* and *Dp* respectively, along the tailings flow.

The increment of the flow mobility caused by the bed material entrainment is clearly shown in the temporal evolution of the flow rate at the control sections CS-2 and CS-3 (**Figure 10**). At CS-2 the increment of the entrainment relation *αE=αD* is mainly appreciated in the increment of the flow rate peak from 7500 m3*=*s for the

*Efficient Simulation Tools (EST) for Sediment-Laden Shallow Flows DOI: http://dx.doi.org/10.5772/intechopen.107345*

#### **Figure 9.**

*Bulk solid concentration in the flow at t* ¼ 25 min *after the dam collapse for (a) αE=αD* ¼ 1%*, (b) αE=αD* ¼ 5%*, (c) αE=αD* ¼ 10% *and (d) αE=αD* ¼ 50%*.*

#### **Figure 10.**

*Temporal evolution of the flow rate at (left) the control sections CS-2 and (right) the control section CS-3 as αE=αD is increased.*

fixed-bed simulation to almost 15,000m3*=*s for the movable bed case with *αE=αD* ¼ 50%. Furthermore, at CS-3, the entrainment of bed material does not only increase the flow rate peak respect to the fixed-bed assumption, but also leads to a noticeable decrease in the arrival time of the mud wave from 40 min to 25 min approximately.

Finally, the computational efficiency of the GPU-accelerated kernel is shown in **Table 2** for different GPU devices. The performance obtained with the GPU-based code is compared with the single-core CPU simulation time. Note that the speed-up obtained respect to the single-core performance has increased progressively as the devices have been developed in the last 5 years. Currently, the speed-up obtained with

#### **Table 2.**

*Computational times with GPU-based and CPU-based algorithms.*

the GPU accelerated kernel is about 280 (Nvidia A100), more than 8 times the Nvidia Tesla K40c performance. That means that for achieving the GPU-parallelized code performance with a CPU-based algorithm, a cluster with at least 280 CPU cores is required (see speed-up in **Table 2**).

#### **5. Conclusions**

The numerical modeling of geophysical surface flows involving sediment transport must be addressed using a comprehensive strategy: beginning with the derivation of proper shallow-type mathematical models, following by the development of robust and accurate numerical schemes within the Finite Volume (FV) framework and ending with the implementation of efficient HPC algorithms. This integrated approach is required to release Efficient Simulation Tools (EST) for environmental processes involving sediment transport with realistic temporal and spatial scales.

The derivation of the generalized 2D system of depth-averaged conservation laws for environmental surface flows of water-sediment mixtures over movable bed conditions is a fundamental step to understand the physical consequences of the mathematical shallow-type simplification. From the mathematical modeling approach, the fluidized material in motion is contained in a flow layer consisting of a mixture of water and multiple solid phases. This flow layer usually moves rapidly downstream steep channels and involves complex topography. The liquid-solid material can be exchanged throughout the bottom interface with the underlying static bed layer, hence involving also a transient bottom boundary for the flow layer. These features lead to an increasing complexity for the mathematical simplification of sediment transport surface flows.

The 2D shallow-type system of equations for variable-density multi-grain watersediment flows can be solved using a Finite Volume (FV) methods, supplemented preferable with upwind augmented Riemann solvers which provide a robust and accurate computation of the intercell fluxes even involving highly transient density interfaces and ensure the well-balanced character of the solution in quiescent and steady states. The usage of GPU-accelerated kernels for sequential computation of the FV solution is an efficient and low cost alternative to the traditional multi-CPU strategies. The GPU solution must be designed taking into account the fact that the GPU is an independent device with its own RAM memory. This means that the memory transfer between the conventional RAM memory and the GPU device memory plays a key role in the performance of GPU-accelerated software.

#### *Efficient Simulation Tools (EST) for Sediment-Laden Shallow Flows DOI: http://dx.doi.org/10.5772/intechopen.107345*

GPU-accelerated codes are mandatory for large-scale sediment-laden flow models since they allows the calibration of the multiple processes involved in the movement of the fluidized water-sediment material without spending weeks or even months waiting results. Probably, the most challenging and unknown process involved in sediment-laden surface flows is the development of pore-fluid pressure, hence its estimation requires careful calibration procedure. This pore pressure affects the frictional shear stress between solid particles, reducing the basal resistance and affecting the flow mobility. Furthermore, special attention has to be paid to the calibration of the entrainment/deposition parameters. The entrainment of material from the underlying movable bed to the mixture layer increases the mass and momentum of the fluidized material in movement, enlarging the mobility of the flow and reducing the arrival time of the mud waves.

#### **Acknowledgements**

This work is part of the research project PGC2018-094341-B-I00 Fast and accurate tools for the simulation and control of environmental flows (FACTOOL) funded by the Ministry of Science and Innovation and FEDER.

### **Author details**

Sergio Martínez-Aranda\*† and Pilar García-Navarro† Fluid Dynamic Technologies (i3A), University of Zaragoza, Zaragoza, Spain

\*Address all correspondence to<EMAIL_ADDRESS>
† These authors contributed equally.

© 2022 The Author(s). Licensee IntechOpen. This chapter is distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

### **References**

[1] Iverson RM. The physics of debris flows. Reviews of Geophysics. 1997; **35**(3):245-296

[2] Pierson TC. Hyperconcentrated flow— Transitional process between water flow and debris flow. In: Debris-flow Hazards and Related Phenomena. Berlin, Germany: Springer Berlin Heidelberg; 2005

[3] Hungr O, Evans SG, Bovis MJ, Hutchinson JN. A review of the classification of landslides of the flow type. Environmental and Engineering Geoscience. 2001;**7**(3):221-238

[4] Rickenmann D, Weber D, Stepanov B. Erosion by debris flows in field and laboratory experiments. In: 3rd International Conference on Debris-Flow Hazards: Mitigation, Mechanics, Prediction and Assessment. Rotterdam: Millpress; 2003. pp. 883-894

[5] Iverson RM, Reid ME, Logan M, LaHusen RG, Godt JW, Griswold JP. Positive feedback and momentum growth during debris-flow entrainment of wet bed sediment. Nature Geoscience. 2011;**4**:116-121

[6] Berger C, McArdell BW, Schlunegger F. Direct measurement of channel erosion by debris flows, illgraben, Switzerland. Journal of Geophysical Research: Earth Surface. 2011;**116**(F1):93-104

[7] McCoy SW, Kean JW, Coe JA, Tucker GE, Staley DM, Wasklewicz TA. Sediment entrainment by debris flows: In situ measurements from the headwaters of a steep catchment. Journal of Geophysical Research: Earth Surface. 2012;**117**:F03016

[8] Wu W. Computational River Dynamics. NetLibrary, Inc. London: CRC Press; 2007

[9] Murillo J, García-Navarro P. Wave riemann description of friction terms in unsteady shallow flows: Application to water and mud/debris floods. Journal of Computational Physics. 2012;**231**: 1963-2001

[10] Juez C, Murillo J, García-Navarro P. 2d simulation of granular flow over irregular steep slopes using global and local coordinates. Journal of Computational Physics. 2013;**255**: 166-204

[11] Gualtieri C, Ianniruberto M, Filizola N. On the mixing of rivers with a difference in density: The case of the negro/solimões confluence, Brazil. Journal of Hydrology. 2019;**578**: 124029

[12] Cunge JA, Holly FM, Verwey A. Practical aspects of computational river hydraulics. In: Monographs and Surveys in Water Resources Engineering. London: Pitman Advanced Publishing Program; 1980

[13] Dewals B, Rulot F, Erpicum S, Archambeau P, Pirotton M. Sediment Transport, Chapter Advanced Topics in Sediment Transport Modelling: Nonalluvial Beds and Hyperconcentrated Flows. Rijeka, Croatia: IntechOpen; 2011

[14] Lacasta A, Juez C, Murillo J, García-Navarro P. An efficient solution for hazardous geophysical flows simulation using GPUs. Computers & Geosciences. 2015;**78**:63-72

[15] Ming X, Liang Q, Xia X, Li D, Fowler HJ. Real-time flood forecasting based on a high-performance 2d hydrodynamic model and numerical weather predictions. Water Resources Research. 2020;**56**:e2019WR025583

#### *Efficient Simulation Tools (EST) for Sediment-Laden Shallow Flows DOI: http://dx.doi.org/10.5772/intechopen.107345*

[16] Juez C, Lacasta A, Murillo J, García-Navarro P. An efficient GPU implementation for a faster simulation of unsteady bed-load transport. Journal of Hydraulic Research. 2016;**54**(3):275-288

[17] de la Asuncion M, Castro MJ. Simulation of tsunamis generated by landslides using adaptive mesh refinement on gpu. Journal of Computational Physics. 2017;**345**:91-110

[18] Pastor M, Blanc T, Haddad B, Drempetic V, Sanchez-Morles M, Dutto P, et al. Depth averaged models for fast landslide propagation: Mathematical, rheological and numerical aspects. Archives of Computational Methods in Engineering. 2015;**22**: 67-104

[19] Jakob M, Hungr O. Debris-Flow Hazards and Related Phenomena. Springer Praxis Books. Springer Berlin Heidelberg; 2005

[20] Iverson RM, Vallance JW. New views of granular mass flows. Geology. 2001;**29**(2):115-118

[21] Berti M, Simoni A. Experimental evidences and numerical modeling of debris flow initiated by channel runoff. Landslides. 2005

[22] McArdell BW, Bartelt P, Kowalski J. Field observations of basal forces and fluid pore pressure in a debris flow. Geophysical Research Letters. 2007; **34**(7):L07406

[23] Richardson JF, Zaki WN. Sedimentation and fluidisation: Part I. Chemical Engineering Research and Design. 1997;**75**:82-100

[24] Godlewski E, Raviart P-A. Numerical Approximation of Hyperbolic Systems of Conservation Laws. New York: Springer-Verlag; 1996

[25] Murillo J, Navas-Montilla A. A comprehensive explanation and exercise of the source terms in hyperbolic systems using Roe type solutions. Application to the 1D-2D shallow water equations. Advances in Water Resources. 2016;**98**:70-96

[26] Murillo J, García-Navarro P. Energy balance numerical schemes for shallow water equations with discontinuous topography. Journal of Computational Physics. 2012;**236**:119-142

[27] Martínez-Aranda S, Murillo J, Morales-Hernández M, García-Navarro P. Novel discretization strategies for the 2D non-Newtonian resistance term in geophysical shallow flows. Engineering Geology. 2022;**302**:106625

[28] Martínez-Aranda S, Murillo J, García-Navarro P. A robust twodimensional model for highly sedimentladen unsteady flows of variable density over movable beds. Journal of Hydroinformatics. 2020;**22**(5):1138-1160

[29] Martínez-Aranda S, Murillo J, García-Navarro P. A GPU-accelerated efficient simulation tool (EST) for 2D variable-density mud/debris flows over non-uniform erodible beds. Engineering Geology. 2021;**296**:106462

[30] Robertson PK, de Melo L, Williams DJ, Ward Wilson G. Report of the Expert Panel on the Technical Causes of the Failure of Feijão Dam I. Technical Report. Brasil: Vale S.A.; 2019

[31] Egiazaroff IV. Calculation of nonuniform sediment concentrations. Proceedings of the American Society of Civil Engineers. 1965;**91**:225-247

### *Edited by Davide Pasquali*

This book presents the latest developments and research in sediment transport. It includes five chapters that address sediment transport estimation, experimental investigation of the physical properties of porous media, and numerical approaches to sediment transport modelling.

Published in London, UK © 2022 IntechOpen © pictureguy32 / Dollarphotoclub

Modeling of Sediment Transport

Modeling of Sediment

Transport

*Edited by Davide Pasquali*


*3.1. Study Sites*

This paper is based on a larger DPhil study that examined the dynamics of land tenure and sustainable land management in Ghana [22]. The analysis presented is based on data collated from two study sites, Kakum and Ankasa Conservation Areas situated in the central and western regions of southern Ghana, respectively (see Figure 1). These landscapes are dominated by the permanently protected Kakum National Park and Ankasa National Park and surrounding communities, spanning a total area of 360 km2 and 509 km2, respectively. The research was conducted in 19 fringe communities randomly selected in Kakum and Ankasa between December 2013 and September 2015 (Figure 1). Aside from having similar tenurial and social contexts affording comparability, both study sites were selected because of their representativeness in illustrating and exploring the prevailing tenurial situation in Ghana's high-forest zone, which spans the southern third of Ghana. In this context, most of the land in the off-reserve areas is stool land wholly owned and managed by the traditional authorities, with pockets of farmlands owned on a leasehold basis (spanning 2 to 50 years) also present [8,22]. Cocoa, oil palm and food crop farming are the dominant economic activities undertaken by households on relatively smaller plots (<5 ha) in the studied communities, working the land under diverse tenurial arrangements, ranging from customary freehold to customary licenses. The tenurial and ethnic diversity in the studied communities afforded a unique opportunity to explore how differences in socio-cultural dynamics (re-)shape tenurial outcomes and conditions of land rights in the study areas.

**Figure 1.** Overview of the two study sites in Ghana. The green area in the map of Ghana illustrates the high-forest zone, with left and right dark arrows pointing to the approximate locations of Ankasa and Kakum conservation areas in the south-western and central regions, respectively. Source: Based on Asaaga and Hirons [8].