Modelling of Sediment Exchange between Suspended-Load and Bed Material in the Middle and Lower Yellow River, China

Abstract

The focus of this paper is on studying novel approaches to estimate sediment exchange between suspended-load and bed material in an unsteady sediment-laden flow with fine-grained sand. The erosion-deposition characteristics of the channel have close relation with the variation of size compositions of both suspended-load and bed material. These aims are addressed by deducing the sediment exchange equations from the mass conservation perspective and establishing a river-sediment mathematical model based on the theory. The model is applied in the middle and lower Yellow River, China, and calibrated and verified under both deposition and erosion conditions using a generalized channel and a large quantity of measured data in the Yellow River basin. The results indicate that the grading curves of suspended-load and bed material calculated by the mathematical model are close to those of the measured data. The temporal and spatial variations in the mean sizes of suspended-load and bed material, flow rate, sediment concentration and erosion or deposition volume estimates during the entire flood process can be accurately predicted. The model performance is considered acceptable for determining the sediment exchange process and the change in channel morphology for unsteady sediment-laden flow.

1. Introduction

Rivers are important sources of water on the earth which are continuously changing, and sediment exchange is one of the most powerful agents in the river environment [1,2]. Sediment exchange transfigures the river morphology by long-term degradation and aggradation of channel beds via erosion and deposition. Sediments in a natural river, including suspended-load, bed-load, and bed material, are usually non-uniform and their size compositions strongly influence the flow capacity and consequent riverbed evolution [3,4,5]. In most alluvial rivers with fine-grained sand in plains (e.g., the middle and lower Yellow River), the proportion of bed-load is low and the sediment exchange between suspended-load and bed material dominates channel morphology change [6,7]. The variations in size compositions of suspended-load and bed material are closely related to erosion-deposition characteristics of the channel [8,9]. The erosion-deposition process is rather complex, which is affected by the interactions between turbulence structures and initiation of particles movement [10]. During an erosion or deposition process, particle coarsening or fining may occur in the same reach simultaneously, which results in the size gradation change of suspended-load and bed material even in the quasi-equilibrium condition, i.e., the amounts of erosion and deposition in the reach can balance each other [11]. Such processes can have a direct effect on the evolution of river forms during flooding. A change in the morphology of a river can threaten channel stability, which can create local scours around hydraulic installations [12]. Evaluating sediment exchange between suspended-load and bed material, therefore, is critical for the long-term prediction of channel morphology change for alluvial rivers with fine-grained sand [13,14,15].

Sediment exchange in alluvial rivers has gained much attention by researchers for years and has been extensively studied using various approaches [16,17]. Among other research methods, laboratorial experimenting is traditional and effective [18]. Jr and Mehta [19] carried out an experiment in a counter-rotating annular flume to explore the exchange between suspended sediment and cohesive bed material. Wang et al. [20] studied the relationship between the bed shear stress and suspended sediment concentration in an annular flume. Zanden et al. [17] presented novel insights into suspended sediment concentrations and fluxes under a large-scale laboratory breaker bar. Qin et al. [14], Yossef and Vriend [21] explored various mechanisms of sediment exchange between the main channel and the groyne fields in a mobile-bed laboratory flume.

Many statistical models for sediment exchange were developed in recent years. Aalto et al. [22] investigated the processes of sediment exchange and channel migration rates along the Strickland River through sampling and Geographic Information System (GIS) analysis. Sinnakaudan et al. [23] focused on sediment transport in high gradient rivers and developed a total-load formula by using the multiple linear regression model [24]. Tayfur et al. [25] and Pektaş [26] used principal component analysis and cluster analysis to identify some significant non-dimensional parameters of sediment transport. Moreover, machine learning based models such as neural networks, fuzzy logic, and support vector machines, were also developed to estimate sediment exchange in recent studies [27,28]. All those machine learning models have a black box nature with limited information inside the models, therefore, they can hardly be expressed as explicit formulas.

Undoubtedly, numerical studies have significantly advanced the understanding of sediment exchange processes and the ability to predict sediment transport rates in alluvial rivers [29]. There has been a plethora of such mathematical models over the last few decades. Many 1D models were developed to simulate the size variations of bed sediment [30,31,32], suspended-load [33], and the hyper-concentration flow [34,35] in alluvial channels. Rudorff et al. [36] predicted the variation in channel-floodplain suspended sediment exchange along a 140 km reach of the lower Amazon River for two decades (1995–2014) with a 2D hydraulic model. Jing et al. [29] proposed a 2D depth averaged Renormalization Group (RNG) k–ε sediment model to simulate the long- and short-term sediment transport and bed deformation in rivers with continuous bends. Aziaian et al. [37], Qamar and Baig [38] applied the CCHE2D model to simulate the non-uniform sediment transport and channel migration in unsteady flows. Yang et al. [39] refined a 2D depth-integrated numerical model based on the Depth Integrated Velocity and Solute Transport (DIVAST) model to simulate bed evolution in the Yellow River Delta, China, from 1992 to 1995. Li and Duffy [40] put forward a 2D high-order model by fully coupling equations of shallow water flow, non-equilibrium sediment transport, and morphological bed evolution in rivers and floodplains. Recently, several 3D models were also developed to simulate sediment movement. Delft3D is capable of modelling 3D hydrodynamics and sediment transport using a finite difference method [11]. Van Rijn et al. [41] took advantage of DELFT3D to simulate the bed evolution of an artificial sand ridge near Hoek van Holland, located close to the harbor of Rotterdam in the Netherlands. Chung and Eppel [42] examining the sensitivity of sediment transport and bed morphology to the variations of bed slopes, grain sizes, and the non-hydrostatic pressure term based on a 3D hydrodynamic and sediment transport model. Franz et al. [43] developed a comprehensive 3D morphological model to predict suspended sediment transport and bed evolution for the Tagus Estuary, Portugal, based on the local hydrodynamic conditions.

To simulate the sediment exchange, most of the above numerical models divide sediment particles into several size groups, treat each group with various methods, and assign some exchanging rules between groups [44]. For instance, Jha and Bombardelli [45] addressed the transport of multi-disperse suspended sediment mixtures in open channels, via the use of the two-fluid model. They confirmed the volumetric concentration of the sediment corresponding to several particle size classes. Bui and Rutschmann [46] developed a sediment transport module for fractional sediment transport using a multiple layer model with the size-fraction method that takes into account the influence of size distribution of the grains on bed-surface on the flow field and bed evolution. Caliskan and Fuhrman [47] coupled a 1D vertical turbulence-closure flow model with sediment transport capabilities to incorporate graded sediment mixtures by treating the bed and suspended load individually for each grain fraction. Franz et al. [43] utilized a morphological model to estimate the non-uniform sediment exchange in natural systems composed of multiple sediment classes. Considering the interactive influences between particles of different sizes have not been well understood, these types of numerical models are always confronted with more complicated computational difficulties and bring considerable numerical inaccuracy [45].

This paper aims to gain new insights into the principles of sediment exchange during riverbed evolution of the alluvial rivers with fine-grained sand, and proposes a novel river-sediment mathematical model to calculate the grain size, composition and exchange of both suspended-load and bed material with no need to divide the sediment into groups based on particle force analysis and mass conservation. The model is applied in the context of the middle and lower Yellow River, China. Although the annual runoff of the Yellow River basin is only approximately 2% of China’s total, it directly supports 12% of the national population and feeds 15% of the irrigation area [48]. Catastrophic floods and droughts appeared frequently in the Yellow River basin over thousands of years in Chinese history, so a large quantity of flow measurements and sediment samples have been collected along the river and stored in a database. This database provides a useful means for validating the model performance. In addition, in the middle and lower Yellow River, bed-load is rare, and the sediment exchange between suspended-load and bed material dominates changes in bed morphology [6,7]. Given its global importance as a large river, the Yellow River is an important case example to evaluate the strengths and potential weaknesses of the theoretical approach for sediment exchange between suspended-load and bed material in alluvial rivers with fine-grained sand presented in this paper.

In this paper, we first present a theoretical framework and establish a river-sediment mathematical model for unsteady flow. Then, a generalized channel is set up to qualitatively analyze and calibrate the model for both deposition and erosion conditions. Finally, the proposed model is applied in the eroded and silted river reaches of the middle and lower Yellow River to investigate its performance. Good agreement between the numerically simulated results and the field measurements is obtained, indicating that the newly-developed model is capable of predicting the sediment gradation and exchange in the middle and lower Yellow River with reasonable accuracy.

2. Study Site and Data

The Yellow River, also called the Huanghe in Chinese, is the second longest river in China. This river is located between 96–119° E longitude and 32–42° N latitude, with a drainage area of c. 753,000 km² and a mean stem length of 5464 km (Figure 1). The river originates from the world’s largest plateau (the Tibetan Plateau), wanders through northern semiarid regions, crosses the loess plateau, passes through the eastern plain, and finally discharges into the world’s largest ocean (the Pacific Ocean) [49]. As seen in Figure 1, the Yellow River basin is usually divided by its physical characteristics into three water source areas: the upper (above Hekou), middle (between Hekou and Huayuankou) and lower (below Huayuankou) reaches [50].

The Yellow River is noted for its relatively low water discharge compared with its very high sediment load: the river has a mean annual water discharge of only approximately 0.7% of the Amazon (the largest river in the world) and 4.5% of the Yangtze (the largest river in China) but an annual sediment load almost equal to that of the Amazon and more than twice that of the Yangtze [51]. The sediment load from the Yellow River basin is estimated to contribute approximately 6% of the total global river load to the oceans [50]. Due to its high sediment concentration and rapid sedimentation rate, the average sediment deposition in the lower reach was 1.58 × 10⁸ t/a (t represents unit ton and a represents age) from 1950 to 2004 [52].

This study focuses on the middle and lower Yellow River as a typical alluvial reach with fine-grained sand (Figure 1). The flow and sediment data collected at 10 hydrological stations along the reach, e.g., Sanmenxia, Tiexie, Huayuankou, Jiahetan, Sunkou, and Lijin, were obtained from the Yellow River Conservancy Commission (YRCC). The physical characteristics of the middle and lower reaches in the Yellow River basin are presented in

Table 1. Physical Characteristics of the Middle and Lower Reaches in the Yellow River Basin [51].

Water Source Area	Length (km)	Drainage Area (km2)	Altitude Range (m)	Annual Mean Precipitation (mm)	Annual Mean Temperature Range (°C)
Middle Reach	1206	343,751	106–1000	530	8–14
Lower Reach	787	22,726	0–106	670	12–14

.

To investigate the applicability of the theoretical approach of sediment exchange and the river-sediment mathematical model presented in this paper, the modelling calculation was conducted using the data of the flood season from the middle and lower Yellow River in 1977 and 1999 to 2002 [48].

Notably, 1977 was a typical year of high sediment load. At Sanmenxia station, the annual water volume was only 30.7 billion m³, and the annual sediment load was 2.08 billion t. In July 1977, the maximum flow was 8100 m³/s, the maximum measured sediment concentration was 535 kg/m³, and the composition of the suspended-load was coarse at Huayuankou station. The hyper-concentrated flood in the lower Yellow River caused channel erosion and beach deposition. Because the deposition volume of beach areas was far greater than the erosion volume of the main channel, the entire channel experienced substantial sediment deposition. During the process of sediment deposition along the channel, the composition of the suspended-load became increasingly fine grained.

After the flood season in 1999, the river reach from Tiexie to Lijin exhibited bed material gradation at 31 measurement sections. The initial composition of the bed material in the lower section was determined according to the distance between sections using an interpolation method. Based on the measured data over three years, the entire channel of the middle and lower Yellow River experienced erosion. Before the flood season in 2002, bed material gradation was observed in the partial middle and lower Yellow River, and the composition of bed material exhibited obvious coarsening compared with that in 1999.

3. Methodology

3.1. Mathematical Model

Based on force analyses of sediment particles under turbulent flow conditions, the variations in gradation of suspended sediment and bed material are derived, as seen in Appendix A1 in detail. According to the principle of mass conservation, unsteady formulas of suspended sediment and bed material exchange are established for the alluvial rivers, as seen in Appendix A2 in detail. The variations of the key parameter, i.e., the mean size of the entrained particles (d_c), the average settling velocity of non-uniform sediment in muddy water (ω_s) and other factors describing the erosion and deposition processes are investigated as seen in Appendix A3. According to these theories, we establish both 1D and 2D river-sediment mathematical models for unsteady flow in this paper. The basic equations of the models are as follows.

1D river-sediment mathematical model:

$$\begin{array}{l}\mathrm{Continuity} \mathrm{equation} \mathrm{of} \mathrm{flow} \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\frac{\partial A_i}{\partial t}+\frac{\partial Q_i}{\partial x}-q_{Li}=0\\ \mathrm{Dynamic} \mathrm{equation} \mathrm{of} \mathrm{flow}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\frac{\partial Q_i}{\partial t}+\frac{\partial }{\partial x}({\alpha }_{1i}\frac{Q_i^2}{A_i})-{\alpha }_{2i}\frac{Q_i}{A_i}{q}_{Li}+g{A}_i(\frac{\partial Z_{wi}}{\partial x}+\frac{Q_i^2}{K_i^2})=0\\ \mathrm{Continuity} \mathrm{equation} \mathrm{of} \mathrm{sediment} \begin{array}{c}\frac{\partial }{\partial t}(A_i{S}_i)+\frac{\partial (A_i{V}_i{S}_i)}{\partial x}+{\displaystyle \sum _{j=1}^m{K}_{1ij}{\alpha }_{*ij}{f}_{1ij}{b}_{ij}{\omega }_{sij}\left(f_{1ij}{S}_{ij}-S_{*ij}\right)}\\ -S_{Li}{q}_{Li}=0\end{array}\\ \mathrm{Riverbed} \mathrm{deformation} \mathrm{equation}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\frac{\partial Z_{bij}}{\partial t}-\frac{K_{1ij}{\alpha }_{*ij}}{{\rho }^{\prime }}{\omega }_{sij}(f_{1ij}{S}_{ij}-S_{*ij})=0\\ \mathrm{Sediment} \mathrm{exchange} \mathrm{equations}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\begin{array}{c}\begin{array}{l}\frac{\partial (QS{d}_{cp})}{\partial x}+\frac{\partial (AS{d}_{cp})}{\partial t}+{\rho }^{\prime }\frac{\partial (A_0{d}_c)}{\partial t}-q_L{S}_L{d}_L=0\\ \frac{\partial D_{cp}}{\partial t}=\frac{d_c-D_{cp}}{H_c}\frac{\partial Z}{\partial t}\end{array}\end{array}\end{array}\},$$

(1)

2D river-sediment mathematical model:

$$\begin{array}{l}\mathrm{Continuity} \mathrm{equation} \mathrm{of} \mathrm{flow} \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\frac{\partial Z_w}{\partial t}+\frac{\partial }{\partial x}(hu)+\frac{\partial }{\partial y}(hv)=0\\ \mathrm{Momentum} \mathrm{equation} \mathrm{of} \mathrm{flow}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em} \begin{array}{c}\frac{\partial u}{\partial t}+u\frac{\partial u}{\partial x}+v\frac{\partial u}{\partial y}+g\frac{\partial Z_w}{\partial x}+g\frac{\sqrt{u^2+v^2}}{C^{2}h}u=\nu \left(\frac{{\partial }^{2}u}{\partial x^2}+\frac{{\partial }^{2}u}{\partial y^2}\right)\\ \frac{\partial v}{\partial t}+u\frac{\partial v}{\partial x}+v\frac{\partial v}{\partial y}+g\frac{\partial Z_w}{\partial y}+g\frac{\sqrt{u^2+v^2}}{C^{2}h}v=\nu \left(\frac{{\partial }^{2}v}{\partial x^2}+\frac{{\partial }^{2}v}{\partial y^2}\right)\end{array}\\ \mathrm{Continuity} \mathrm{equation} \mathrm{of} \mathrm{sediment} \frac{\partial (hS)}{\partial t}+\nabla \times (hS\overrightarrow{V})=\epsilon {\nabla }^2(hS)+{\alpha }_{\ast }{\omega }_s{K}_1(S_{\ast }-f_{1}S)\\ \mathrm{Riverbed} \mathrm{deformation} \mathrm{equation}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\frac{\partial Z_{bij}}{\partial t}-\frac{K_{1ij}{\alpha }_{*ij}}{{\rho }^{\prime }}{\omega }_{sij}(f_{1ij}{S}_{ij}-S_{*ij})=0\\ \mathrm{Sediment} \mathrm{exchange} \mathrm{equations}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\begin{array}{c}\begin{array}{l}\begin{array}{c}uhS\frac{\partial d_{cp}}{\partial x}+vhS\frac{\partial d_{cp}}{\partial y}-h{\epsilon }_{sx}\frac{\partial S}{\partial x}\frac{\partial d_{cp}}{\partial x}-\frac{\partial }{\partial x}\left(hS{\epsilon }_{sx}\frac{\partial d_{cp}}{\partial x}\right)\\ -h{\epsilon }_{sy}\frac{\partial S}{\partial y}\frac{\partial d_{cp}}{\partial y}-\frac{\partial }{\partial y}\left(hS{\epsilon }_{sy}\frac{\partial d_{cp}}{\partial y}\right)+hS\frac{\partial d_{cp}}{\partial t}\\ +{\rho }^{\prime }\frac{\partial \left(Z_S{d}_c\right)}{\partial t}-{\rho }^{\prime }{d}_{cp}\frac{\partial Z_S}{\partial t}=0\end{array}\\ \frac{\partial D_{cp}}{\partial t}=\frac{d_c-D_{cp}}{H_c}\frac{\partial Z}{\partial t}\end{array}\end{array}\end{array}\},$$

(2)

where subscripts “i” and “j” represent each cross section and sub-cross section, respectively, m is the number of sub-cross sections, Z_w is the water level, K is the discharge modulus, α₁ and α₂ are the momentum correction coefficients of flow and lateral inflow, respectively, α_* is the distribution coefficient of sediment concentration in an equilibrium state, Z_b is the average bed elevation of each cross section, C is Chezy coefficient, S_* is the vertical averaged sediment carrying capacity of the flow, and d_cp and D_cp are the mean sizes of total suspended sediment and bed material, respectively. Full details of the other definitions of the model are presented in the Appendix A. The primary parameters during the sediment exchange process in above models can be ascertained as follows.

$${\omega }_s={\omega }_0(1-1.25S_V){[1-\frac{S_V}{2.25\sqrt{d_{50}}}]}^{3.5},$$

(3)

$$S_{*}=2.5{\left[\frac{\left(0.0022+S_V\right)V^3}{\kappa \frac{{\gamma }_s-{\gamma }_m}{{\gamma }_m}gh{\omega }_s}\ln \left(\frac{h}{6D_{50}}\right)\right]}^{0.62},$$

(4)

$$d_c=\{\begin{array}{l}{\left(\frac{\Delta Z{\rho }^{\prime }}{Sh}\right)}^{\frac{d_{cp}}{d_{ph}-d_{cp}}}{d}_{cp}+\left[1-{\left(\frac{\Delta Z{\rho }^{\prime }}{Sh}\right)}^{\frac{d_{cp}}{d_{ph}-d_{cp}}}\right]d_{ph}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\mathrm{Deposition}\\ {\left(\frac{\left|\Delta Z\right|}{{H^{\prime }}_c}\right)}^{\frac{D^{″}-d_f}{d_f-d_{ph}}}{d}_f+\left[1-{\left(\frac{\left|\Delta Z\right|}{{H^{\prime }}_c}\right)}^{\frac{D^{″}-d_f}{d_f-d_{ph}}}\right]d_{ph}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\mathrm{Erosion}\end{array},$$

(5)

where ${\omega }_0=2.6{(\frac{d_{cp}}{d_{50}})}^{0.3}{\omega }_p\times e^{-635d_{cp}^{0.7}}$, which is the average settling velocity of non-uniform sediment in clear water; $d_{ph}=N_d{d}_{cp}\exp \left\{8.21\left[f\left(n\right)-1\right]\frac{{\omega }_s}{\kappa u_{\ast }}\right\}$, which is the mean size of bed material load when channel morphology change reaches dynamic equilibrium; $d_{f\max }=M{\left(\frac{\gamma }{{\gamma }_s-\gamma }\frac{1}{g}\right)}^{2/3}{\left(\frac{D_{cp}}{h}\right)}^{1/6}{\nu }^{1/3}V$, which is the maximum particle diameter under critical conditions for sediment suspension; $D^{″}$ is the particle diameter of relatively coarse bed material; and M is a coefficient.

3.2. Modelling Procedure and Scenario Setting

When solving the mathematical models, the uncoupled solution is adopted. Firstly, solve the continuity and dynamic equations of flow and obtain the hydraulic elements. Then, solve the continuity equation of sediment and riverbed deformation equation, and ascertain the erosion and deposition process of the riverbed. Finally, calculate d_cp and D_cp. The solving process proceeds alternately and circularly till the simulation period is over.

Two simulation stages are designed in this study. The first stage is to calibrate the presented models. In the second stage, abundant simulations of a hyper-concentrated flood of the lower Yellow River are conducted to evaluate the applicability of the mathematical models of sediment exchange for both the deposition and erosion conditions.

3.3. Model Calibration and Verification

To preliminarily calibrate the model presented in this paper, we carry out qualitative analysis and calculation of a generalized channel. According to the on-site sediment composition and river bed gradient in the middle and lower Yellow River [48], the generalized channel is set as 100 km long and 200 m wide with a vertical gradient of 1$\raisebox{1ex}{$0$}\!\left/ \!\raisebox{-1ex}{$000$}\right.$. The initial mean and median sizes of the bed material (D_cp and D₅₀) are set at 0.125 mm and 0.1 mm, respectively. The roughness coefficient of the riverbed is set at 0.012. The steady flow rate is set at 530 m³/s, and the initial water depth is 2 m. For the deposition conditions, the sediment concentration S of the inlet is set at 100 kg/m³, and the mean and median sizes of the suspended-load (d_cp and d₅₀) are set at 0.04 mm and 0.032 mm, respectively. For the erosion conditions, the sediment concentration S of the inlet is set at 1 kg/m³, and the mean and median sizes of the suspended-load (d_cp and d₅₀) are set at 0.015 mm and 0.012 mm, respectively. The temporal and spatial variations in the mean and median sizes of suspended-load and bed material of the generalized channel are obtained by the model simulation.

Moreover, we have derived sediment gradation equations in light of force analyses of sediment particles (Appendix A1). The sediment gradation can be calculated based on the simulation results to verify the mathematical models. The sediment gradation equations are as follows:

$$\{\begin{array}{l}P(d_i)=2\mathsf{\Phi }\left[0.675{(\frac{d_i}{d_{50}})}^n\right]-1\\ n=0.42{\left[\mathrm{tg}\left(1.49\frac{d_{50}}{d_{cp}}\right)\right]}^{0.61}+0.143\end{array},$$

(6)

$$\{\begin{array}{l}P(D_i)=2\mathsf{\Phi }\left[0.675{(\frac{D_i}{D_{50}})}^n\right]-1\\ n=0.42{\left[\mathrm{tg}\left(1.49\frac{D_{50}}{D_{cp}}\right)\right]}^{0.61}+0.143\end{array},$$

(7)

where $\mathsf{\Phi }$ is the normal distribution function.

To evaluate the performance of the newly developed model, the calculated and measured values are compared using the Nash-Sutcliffe efficiency coefficient (NSE) and Pearson correlation coefficient r defined as follows:

$$NSE=1-\frac{{\displaystyle \sum _{i-1}^n{\left(X_{m,i}-X_{c,i}\right)}^2}}{{\displaystyle \sum _{i-1}^n{\left(X_{m,i}-{\overline{X}}_m\right)}^2}},$$

(8)

$$r=\frac{{\displaystyle \sum _{i-1}^n\left(X_{m,i}-{\overline{X}}_m\right)\left(X_{c,i}-{\overline{X}}_c\right)}}{\sqrt{{\displaystyle \sum _{i-1}^n{\left(X_{m,i}-{\overline{X}}_m\right)}^2}{\displaystyle \sum _{i-1}^n{\left(X_{c,i}-{\overline{X}}_c\right)}^2}}},$$

(9)

where X_m and X_c are the measured and calculated values, respectively. NSE and r are commonly used to quantitatively assess the predictive power of hydrological models [53,54]. If NSE > 0.5 and r > 0.6, the model performance is considered to be acceptable.

3.4. Model Settings and Parameters

To investigate the applicability of the theoretical approach of sediment exchange and the river-sediment mathematical model presented in this paper, measured data from a hyper-concentrated flood of the lower Yellow River in July 1977 were used in the analysis for the deposition conditions, and the river reach used for verification extended from Huayuankou to Jiahetan (Figure 1). The initial bed material gradation in each section was determined according to the measured gradation data from Huayuankou to Jiahetan before the flood season based on the interpolation method. The outlet boundary condition was controlled by the stage-discharge relationship of Jiahetan station in 1977.

The temporally variable values of the inlet flow and sediment concentration in the model were determined according to the measured data from Huayuankou station using a temporal interpolation method. Because the quantity of measured gradation data for the suspended-load was small, the modeled results differ from the actual results if the temporal interpolation method is directly applied. As a result, a relation involving d_cp, Q and S was established through regression analysis.

$$d_{cp}=0.00056Q^{0.224}{S}^{0.428}$$

(10)

The value of d_cp of the suspended-load corresponding to the measurement time of the sediment concentration was obtained using Equation (10), and the existing measured value of d_cp was retained. Then, the values of d_cp during the entire computational period were determined using the temporal interpolation method.

To verify the newly proposed approaches for erosion conditions in this paper, the data collected from the middle and lower Yellow River during the period between November 1999 and October 2002 were used in this analysis, and the river reach used for verification extended from Tiexie to Lijin (Figure 1). The measured water and sediment values and the composition of the suspended-load at Xiaolangdi station were set as the inlet boundary conditions. The initial topography was controlled by the measured sections of the partial middle and lower Yellow River after the flood season in 1999. The water level at the end of the model was controlled using the stage-discharge relationship observed at Lijin station in the same year.

4. Results and Discussion

4.1. Results of Model Calibration

4.1.1. Deposition Conditions

Figure 2 shows the temporal variations in $d_{cp}$, $d_{50}$ and $d_{50}/d_{cp}$ with different cross sections along the channel. In the silting process, $d_{cp}$ and $d_{50}$ decrease along the channel, and $d_{50}/d_{cp}$ increases along the channel. The degree of sediment fining of the upper reach (30–60 km) is larger than that of the lower reach (60–90 km) over an equal distance and time frame. During the calculated period, the sediment-carrying capacity of the flow gradually increases with the increase in the vertical gradient as a result of deposition, so the fining effect of the suspended-load has a tendency to decrease.

Figure 3 presents the spatial variations in $D_{cp}$, $D_{50}$ and $D_{50}/D_{cp}$ along the channel with different deposition times. As shown, the composition of bed material gradually becomes finer in the silting process. The fining process starts from the upper reach to the lower reach. In the early stage, the bed material becomes more non-uniform as abundant suspended-load silts on the bed ($D_{50}/D_{cp}$ decreases in the same cross section). With increasing deposition, the bed material gradation becomes relatively uniform ($D_{50}/D_{cp}$ increases with a calculation time of 100 days).

4.1.2. Erosion Conditions

As shown in Figure 4, $d_{cp}$ and $d_{50}$ increase along the channel during scouring, and $d_{50}/d_{cp}$ decreases along the channel, which means the composition of the suspended-load becomes increasingly non-uniform. The degree of sediment coarsening of the upper reach is larger than that of the lower reach over the same time. The sediment-carrying capacity of flow gradually decreases with the decreasing vertical gradient of the riverbed during the erosion process, and the coarsening effect of the suspended-load tends to decrease.

As shown in Figure 5, bed scouring starts from the upper reach to the lower reach with the increase in $D_{cp}$ and $D_{50}$ during scouring, and the coarsening effect on the bed material in the upper reach is larger than that in the lower reach. In the early stage of erosion, the bed material becomes increasingly non-uniform as the finer bed material is scoured up to the flow, so $D_{50}/D_{cp}$ decreases in the same cross section. With the development of erosion, more and more fine particles are washed away, so the bed material gradation gradually becomes uniform. As the erosion quantity of the upper reach is larger than that of the lower reach, $D_{50}/D_{cp}$ in some sections of the upper reach begins to increase with the calculation time of 10 days.

In light of the calculation results of the generalized channel under deposition or erosion conditions, the variation in sediment gradation and exchange is consistent with the basic rules of sediment movement. Thus, the model presented in this paper is considered to be preliminarily reliable.

4.2. Application and Verification

4.2.1. Verification of the Sediment Gradation Calculation

We calculated the sediment gradation based on the model and Equations (6) and (7). Figure 6a,b compare the calculated grading curves of the suspended-load and bed material, respectively, with the measurements at three hydrological stations in the lower Yellow River, i.e., Huangyuankou, Lijin, and Gaocun stations. The model can clearly be applied to calculate the suspended-load and bed material gradation of a natural river with a wide range of particle diameters.

4.2.2. Verification of the Sediment Exchange Calculation for Deposition Conditions

Based on the model simulation results, Figure 7 shows the comparative results of the calculated and measured values of the flow rate (NSE = 0.95 and r = 0.95), sediment concentration (NSE = 0.98 and r = 0.98), and d_cp (NSE = 0.89 and r = 0.97). The results of the exchange simulation method of the suspended-load and bed material presented in this paper exhibit good agreement with the measured data for the peak flow, flood transmission time, maximum sediment concentration and its emergence time. The results present the high calculation accuracy of the flow and sediment concentration during the entire process compared to the results of traditional methods. Moreover, the values of and variations in d_cp for the suspended-load are relatively close to the measured values and trends. The established mathematical model can further improve the accuracy and functionality of the simulated evolution of hyper-concentrated floods and variations in riverbed deposition in the lower Yellow River. Thus, the method of simulating exchange processes involving suspended-load and bed material presented in this paper is applicable for deposition conditions in natural rivers.

4.2.3. Verification of the Sediment Exchange Calculation for Erosion Conditions

Figure 8 displays the comparative results of the calculated and measured erosion processes for the entire lower channel (NSE = 0.98 and r = 0.98). Figure 9 displays the comparative results of the mean size of the bed material D_cp for the major channel surface of the middle and lower Yellow River before the flood season in 2002 (NSE = 0.78 and r = 0.93). According to the comparative results, the evolution process and distribution of erosion in the middle and lower Yellow River calculated by the model agree well with the measured data. Additionally, the coarsening of bed material by erosion and exchange between bed material and suspended-load are supported by the results. The computational accuracy of the model satisfies the production requirement. As a result, the method of simulating the exchange process of suspended-load and bed material presented in this paper is applicable under erosion conditions in natural rivers.

5. Conclusions

In this paper, we have systematically presented a novel approach for simulating the gradation and exchange of suspended sediment and bed material for the alluvial rivers with fine-grained sand. Based on force analyses of sediment particles under turbulent flow conditions, the variations in sediment gradation are calculated. According to the principle of mass conservation, 1D, 2D and 3D unsteady formulas of suspended sediment and bed material exchange are established. This approach avoids the adverse factors compared with the calculation methods based on the fractional carrying capacity for graded sediment, and a river-sediment mathematical model is built in light of the theoretical approach. The model is applied to the case of the middle and lower Yellow River, China, which belongs to representative unsteady sediment-laden flow with high concentration of fine-grained sand. Some of the most important river-sediment features during flood events are reproduced, including the sediment particle diameter, flow rate, sediment concentration, erosion and deposition volume, etc.

The results imply that the grading curves of the suspended-load and bed material calculated are relatively close to those of the measured data. The predicted values of the mean sizes of the suspended-load and bed material match well with the associated measured values, and the computational accuracy of the flow, sediment concentration and erosion and deposition volume estimates during the entire flood process is greatly improved. In contrast with the traditional methods, the newly developed approaches could not only provide an effective tool for simulating the sediment exchange between suspended-load and bed material for sediment-laden flow, but also improve the ability of numerical models in estimating sediment concentration and channel morphology change during flood events.