Experimental and numerical study of bed roughness effect on longitudinal dispersion

When the solutes are discharged to the river systems, and as they are transported downstream by the flowing water, they will encounter different hydraulic and geochemical processes, such as mixing, exchange with storage zones and biogeochemical reactions. The effluent cloud initially spreads vertically, transversely, and longitudinally due to turbulent and molecular diffusion and shear dispersion. The dispersion mechanisms in rivers are due to shear velocity and due to storage zone effects (Davis et al. 2000; Atkinson & Davis 2005). The first one is due to non-uniform velocity distribution and turbulent diffusion and it may lead to a Fickian distribution of the solute concentration. The dispersion due to the storage zone denotes tracer exchange between the shear flow and the storage zones, in which, there is no longitudinal transport (Chen et al. 2016). Deng et al. (2001) divided dispersion into three processes of molecular diffusion, turbulent diffusion, and shear dispersion. They pointed out that the shear dispersion is induced by velocity gradients due to the fluid viscosity and the resistance of the fluid along boundaries. Shear dispersion occurs in both laminar and turbulent flows and is dominant compared to the diffusion mechanisms.

In river systems, the river bed usually contains non-uniform geometries (for example, sharp bends and meandering structures, sinuosity, cross-section shape, riparian vegetation, suspended sediment, and water surface waves) (Ng 2000; Marion & Zaramella 2006; Boxall & Guymer 2007; Perucca et al. 2009; Huang & Law 2011). Therefore, the bed roughness may greatly affect key factors contributing to the longitudinal dispersion such as shear stress distribution, velocity distribution, vertical turbulent diffusivity, and transverse mixing. Subsequently, one can expect a different longitudinal dispersion coefficient in rough bed rivers from that in smooth bed streams. Conversely, longitudinal dispersion coefficients (D_L) often show different values for the same river reach and discharge (Seo & Cheong 2001; Deng & Jung 2009; Shen et al. 2010). As D_L cannot be measured directly, and tracer injection tests for indirect estimation of D_L are often time consuming and costly, empirical equations have been widely employed to calculate D_L for natural open channels. These equations have mainly been formulated based on some available hydraulic-geometric parameters, such as the cross-sectional averaged velocity (U), shear velocity (U*), flow depth (d) and channel width (W) (Kashefipour & Falconer 2002; Noori et al. 2011; Zeng & Huai 2014; Balf et al. 2018). Tenebe et al. (2016) reviewed the various models derived as well as methods associated in the collection of tracer concentration data existing in the literature.

Many experimental works have been carried out on the flow dispersion over rough bed surfaces. Elder (1959) and Fischer et al. (1979) studied the longitudinal dispersion coefficients in natural streams. Fischer (1967) reported dispersion coefficients in the range of 0.021–0.047 m²s⁻¹ in small trapezoidal artificial channels with roughened sidewalls. Abd El-Hadi & Davar (1976) evaluated the longitudinal dispersion coefficient (D_L) over the large bed roughness in a laboratory flume. Their results showed that the dimensionless dispersion coefficient, D_L/(hU*) (in which, h is the water depth and U* is the bed shear velocity), is a function of the relative roughness height (k_s/h) and the relative roughness spacing. It was also found that D_L/(hU*) reaches a maximum value at a relative roughness spacing of 0.111. Magazine (Magazine 1983) analyzed the dimensionless dispersion coefficient as a function of two interacting mechanisms of resistance and blockage effects due to roughness elements. They concluded that the dispersion characteristics of rough bed channels are analogous to those of rough side channels. They also proposed the parameter P, incorporating the effect of resistance and blockage effects and found a good correlation between this parameter and available dispersion data. They also developed an empirical model that can be used to predict dispersion coefficients in the straight laboratory and natural rivers. Seo & Cheong (1998) argued that large spacing and height of roughness might lead to resistance and blockage effects at the fluid boundary. Thereby, the gradients in transverse and vertical velocity distribution are increased, which in turn enhances the longitudinal mixing. Runkel (2002) found that the longitudinal dispersion coefficients in open channels (man-made channels) were considerably smaller (about ten times) than those found in small natural streams. This can be attributed to the uniform geometry of the channel. Deng et al. (2002) presented datasets of 70 and 53 dispersion coefficients, respectively, measured in streams and rivers at various scales and ranging from 1.9 to 1,486 m²s⁻¹. In addition, the large variability of longitudinal dispersion in dependency on the experimental scale was demonstrated. Schulz et al. (2012) investigated the effect of bed roughness on longitudinal dispersion in artificial flow channels by tracer experiments with variations in channel bed material. They found that the longitudinal dispersion coefficients ranged from 0.018 m²s⁻¹ in smooth bed channels to 0.209 m²s⁻¹ in coarse gravel bed channels. They concluded that the longitudinal dispersion was linearly related to the average flow velocity. In addition, longitudinal dispersivities (λ = D_L/V) ranged between 0.152 ± 0.017 m in smooth bed channels and 0.584 ± 0.015 m in the coarse gravel bed channel. Grain size and bed roughness were positively correlated to longitudinal dispersion. They also pointed out that this finding contradicted several existing relations between bed roughness and longitudinal dispersion. Bahadur et al. (2021) compared 30 equations for prediction of the longitudinal dispersion coefficient by using hydrologic data from 59 river reaches. Their results showed that sinuosity significantly impacts estimation of the dispersion coefficient. In addition, the computations that include sinuosity improved the performance of the dispersion equation.

Many researchers have numerically resolved the advection–dispersion equation for idealized cases (Guymon 1970; Smith et al. 1973; Ehlig 1977; Lam 1977). Such numerical methods involve two problems. Firstly, an accurate estimate of D_L must be made, and secondly, a suitable mathematical algorithm must be developed. Deng et al. (2004) developed numerical schemes and stability criteria for solution of the one-dimensional fractional advection–dispersion equation (FRADE) derived by revising Fick's first law. Barati Moghaddam et al. (2017) presented the TOASTS (third-order accuracy simulation of transient storage) model to simulate 1-D pollutant transport in rivers with irregular cross-sections under unsteady flow and transient storage zones. The proposed model was verified versus some analytical solutions and a 2-D hydrodynamic model. The results of the TOASTS model, in comparison with two common contaminant transport models, shows better accuracy and numerical stability. Parsaie & Haghiabi (2017) solved fractional advection–dispersion equation (FRADE) and developed a numerical model for the simulation of pollution transmission in rivers with stagnant zones. To this purpose, both terms of FRADE equation (advection and fractional dispersion) were discreted separately and the results of them were connected together with the time-splitting technique. Their results indicated that there is a good agreement between observed data, the analytical solution of MADE model, and results of developed numerical model. Ramezani et al. (2019) presented a numerical modelling-based approach to evaluate the performance of various D_L formulas using the ADE. This approach was tested against the analytical ADE solution and demonstrated using eight well-known D_L formulas and tracer study data for the Chattahoochee River (USA), the Severn (UK) and the Athabasca (Canada). Their results showed that D_L has an important effect on tracer concentrations simulated with the ADE. Comparison between the simulated and measured concentrations confirmed the appropriate performance of Zeng & Huai (2014)'s formula for D_L estimation.

In addition, some researchers utilized data-mining methods for the estimation of the dispersion coefficients. Sulaiman et al. (2018) developed a hybrid soft computing model called deep neural network (DNN) coupled with a genetic algorithm (GA) to predict LDC using historical information obtained from published researches in the literature. Arya Azar et al. (2021) predicted the longitudinal dispersion coefficient in natural streams using LS-SVM and ANFIS. The results showed that machine learning models had a more efficient performance to predict the longitudinal dispersion coefficient than experimental equations. Ghiasi et al. (2021) used a Deep Convolutional Network (DCN), a sub-field of machine learning for predicting the longitudinal dispersion coefficient.

The above-mentioned studies in surface waters showed that various parameters such as bed roughness and river geometry have a significant influence on the longitudinal dispersion. Although hydraulic engineers have been studying the bed roughness effect for decades, still, there are few comprehensive studies on the effect of bed roughness on longitudinal dispersion in open channels and many well-known developed equations for estimation of longitudinal dispersion do not consider the bed roughness effect (Noss & Lorke 2016; Wang & Huai 2016). Therefore, the main research gap and the reason for doing this research is, first, what effect does bed roughness have on longitudinal dispersion, and second, to what extent does bed roughness change this coefficient, and what is the importance of considering it in the estimation equations. Of interest in this study is experimental and numerical modelling of the advection–dispersion equation at the presence of bed roughness. To achieve an overall output in this regard, the tracer experiments were carried out with a variety of gravels, each having the uniform size distribution. By using the measured tracer breakthrough curves and average flow velocity, the mean longitudinal dispersion coefficient and mean longitudinal dispersivity were numerically calculated. Finally, a relationship was developed between relative roughness and dimensionless longitudinal dispersion coefficient. This research will have a significant contribution to determine the influence of bed roughness on the longitudinal dispersion coefficient. Section 2 deals with the methodology of the conducted experiments and developed numerical model. In section 3, the results are presented and analyzed. Using these results, a regression equation was developed to determine the dispersion coefficient of the flow over rough bed channel.

The experimental work was conducted in a circular flume with a mean diameter of 1.6 m, a width of 0.20 m, and a depth of 0.15 m (Figure 1). The flume was placed on a stationary 2 × 2 m² platform, and flow was run using two rotating pedals within the flume (Mahdavi et al. 2013). Circular flumes give the researcher the opportunity to study the contamination transport in these flumes by changing the boundary conditions, instead of using the straight and long flumes.

A series of tracer experiments (48 runs) was conducted using four different uniform bed materials of 1.04, 2.09, 3.01, and 4.24 mm (coarse sand). According to Einstein criteria, k_s = D₆₅ was used for bed roughness heights (Einstein & El-Samni 1949). A bare plexiglass wall was utilized as a reference to evaluate the effect of bed roughness in comparison with the smooth bed channels. The experiments were performed by laying down and sticking each coarse sand with a uniform distribution over the flume bottom to produce a uniform rough bed. Figure 2 shows the pictures of gravel materials used for creating a rough bed. Also, the grain size distribution and analysis of coarse sands for rough beds are presented in Figure 3 and Table 1, respectively. The experiments were performed at different flow velocities of 0.15 to 0.35 ms⁻¹. These values were selected according to values observed in natural streams. Moreover, for each experiment, the flow depth (h) was kept at a constant height of 13 cm above the flume base.

In these experiments, sodium chloride (NaCl) was used as tracer matter, and two conductivity meters were used to measure the time distribution of the sodium chloride concentration (breakthrough curves). After stabilizing average flow velocity and the environmental parameters (such as initial electrical conductivity and temperature), 250 mL of NaCl solution with a specified concentration was injected from the point x = 0 (see Figure 1). The start time was recorded from this moment. The conductivity meters were installed at a distance of x = 1.256 m from the injection point (see Figure 1), recording NaCl concentration at 2 s intervals in the memory. The list of conducted experiments is presented in Table 2. The average flow velocity was estimated by installing the conductivity meters at a specified point of the flume and measuring the breakthrough curves of the tracer. Figure 4 shows a sample of the tracer breakthrough curve. Given that the solute cloud moves at a rate equal to the average flow velocity throughout the flume, the distance between the two peaks indicates the duration of rotation of clouds within the flume. It should be mentioned that breakthrough curves were recorded for at least seven or eight circulations (Figure 4). As can be seen, the baselines in some breakthrough curves increase, which is attributed to the longitudinal dispersion and transient storage, both of which lead to the overlapping of subsequent peaks (Schulz et al. 2012).

NaCl breakthrough curves were recorded until the solute cloud was completely spread throughout the flume, and the solute concentration was constant at different points of the flume (see Figure 4).

A numerical model was developed by solving Equation (1) based on Apollo and Postman's (Appelo & Postma 1993) method. According to this method, the advection and dispersion terms are individually solved at each step. In other words, at a certain time step, the effect of the advection on the cells is first considered. Then, the effect of the dispersion is applied to the previous result. The advantage of this method is that one can easily add other chemical or physical formulas to the model (Appelo & Postma 1993; Mahdavi et al. 2013). In addition, each process can be solved with the most appropriate numerical method.

In Equation (2), t refers to the time step and i denotes the spatial step. For stability of numerical model, the Neumann number (DΔt/(Δx²)) must be less than 0.5. After determining the time step using the flow velocity and the spatial step, the time step is divided into M steps in such a way that the Newman number is less than 0.5, then the dispersion equation will be repeated for M times.

The initial condition is considered as at time t = 0, the concentration across the flume except for the injection site is zero. In addition, the boundary condition is taken as the concentration at x = 0 and x = L are equal at all times.

As mentioned before, the tracer experiments were conducted to determine average flow velocity and longitudinal dispersion coefficients. In total, 48 tracer experiments were performed at two different flow velocities and two different injection concentrations, as well as two initial concentrations. Flow velocities were kept constant during each experiment; however, they varied between 0.2 and 0.35 ms⁻¹.

Figure 6 shows the time distributions NaCl for the experiments with two different average flow velocities on both smooth and rough beds (k_s = 4.24 mm). In these figures, along with experimental data, the results of numerical solution are also plotted. Solid and dashed lines in these figures are related to numerical simulation of the flow over smooth and rough beds, respectively. In addition, dark and gray points are respectively relevant to the experiments on smooth and rough beds (k_s = 4.24 mm). The numerical model was then calibrated to determine the dispersion coefficient. The model was run for different flow velocities, and the best dispersion coefficient was selected after comparing with the experimental data and calculating the errors based on the statistical indices.

As can be seen in these figures, the tracer concentrations in the smooth bed are considerably different from those in the rough bed. It is also obvious that the distance between the peaks decreases, indicating bed roughness is reduced the average flow velocity. On the other hand, the time taken for the solute cloud to spread completely throughout the flume is significantly decreased over the rough bed. In other words, the breakthrough curves in rough beds have a lower peak concentration and longer that for the smooth bed. Another notable change in these diagrams is that the number of circulations in breakthrough curves in the rough bed are significantly reduced compared with the smooth bed. For example, in Figure 7(d), the number of circulations has decreased from 16 circulations in smooth bed to nine circulations on a rough bed. This is evidence of the increasing dispersion rate of the solute on the rough beds.

The results of experiments on smooth and rough beds with different roughness heights are presented in Table 3. As shown in this table, the bed roughness reduced the average flow velocity from 0.189 to 0.141 ms⁻¹ (up to about 35%, for minimum velocity). In addition, longitudinal dispersion coefficients ranged from 0.00925 m²s⁻¹ in smooth bed experiments up to 0.012 m²s⁻¹ in channels with coarse sands, indicating an increase of about 34%. These changes in average flow velocity and the longitudinal dispersion coefficients are greater for the experiments with the maximum flow velocity. Increasing the bed roughness height may increase the mixing of fluids, especially near the bed. The mixing process may reduce the flow velocity near the bed and increase the longitudinal dispersion coefficient. This shows that the presence of the sediment in the river bed can also increase the longitudinal dispersion coefficient. Furthermore, the increase in longitudinal dispersion coefficient in the rough bed experiments shows that the flow turbulence and mixing are enhanced compared to the smooth bed.

As can be seen in Table 3, with increasing relative roughness, the values of longitudinal dispersivity (λ) are increased. The range of λ in smooth bed experiments is 0.037–0.049 m, and for rough beds (with a roughness height of 4.24 mm) the values of λ are in the range of 0.07–0.084 m. The results showed that, with increasing the roughness height from zero to 4.24 mm, the longitudinal dispersivities increased from 0.037 m to 0.077 m, indicating an increase of about 108%. According to the main aim of the present study, the presence of bed roughness materials may significantly affect the longitudinal dispersion in natural streams.

For the sake of comparison of the results of present study with previous studies, the study of Schulz et al. (2012) was selected as the most relevant study. They conducted the experiment in three modes of no gravel (smooth bed), fine gravel (2–8 mm) and coarse gravel (32–64 mm). The range of average flow velocities in their experiment was 0.306–0.366 ms⁻¹. They calculated the dispersion coefficients using two statistical methods. They concluded that the dispersion coefficients in smooth bed, fine gravel and coarse gravel were respectively 0.048, 0.139, 0.209 m²s⁻¹, an increase of about 335%. In other words, the dispersion coefficients in coarse gravel beds were approximately four times the dispersion coefficients in smooth beds.

The calculated longitudinal dispersion coefficients are plotted against relative roughness in Figure 8 for two different flow velocities. In this figure, it is well obvious that, with increasing the relative roughness, the longitudinal dispersion coefficient is increased. It is also observed that with increasing the flow velocity from about 0.2 to about 0.35 ms⁻¹, the longitudinal dispersion coefficient increased from about 0.001 m²s⁻¹ up to about 0.014 m²s⁻¹.

In order to establish a relationship between the longitudinal dispersion coefficients and the roughness height, the non-dimensional parameter of D_L/(Vh) was plotted against the relative roughness (k_s/h) (Figure 9).

In order to verify the results of this equation, three statistical indices of mean absolute error (MAE), sum of the squared errors (SSE), and determination coefficient (R²) were used. The results of the analysis with these statistical indices showed that MAE, SSE, and R² are respectively 47.7%, 0.07, and 0.88. These values represent the great accuracy of this relationship in estimating the dispersion coefficient in the presence of bed roughness.

The purpose of this study was to investigate the relative effect of different bed roughness on longitudinal dispersion. The quantitative estimates of the errors between the observed and calculated values of D_L indicated that the numerical model developed in the present study was reasonably applicable to compute the tracer concentrations on both the smooth and rough beds. The results showed that the tracer concentrations in the smooth bed are considerably different from those in the rough bed. In addition, the time taken for the solute cloud to completely spread throughout the flume is significantly decreased over the rough bed. It was observed that the bed roughness reduced the average flow velocity up to about 35% and increased the longitudinal dispersion coefficients about 34%. These changes are greater for the experiments with the maximum flow velocity. This shows that the presence of sediment in the river bed can increase the longitudinal dispersion coefficient. Accordingly, with increasing the relative roughness, the values of longitudinal dispersivity increase. The results showed that with increasing the roughness height from zero to 4.24 mm, the longitudinal dispersivities are increased by about 108%. Additionally, a relationship was developed using non-dimensional longitudinal dispersion (D_L/Vh) as a function of relative roughness (k_s/h). The results showed the high ability of this relationship to estimate the longitudinal dispersion coefficient in both smooth and rough beds. It can be concluded that the consideration of bed roughness as the driving force of shear dispersion would improve predictive equations of longitudinal dispersion in the rivers. As the bottom of all natural rivers has roughness elements with different sizes, the results of this study will definitely be useful in estimating the longitudinal dispersion coefficient in natural rivers and quantifying the effect of roughness in the longitudinal dispersion coefficient equations.

One of the limitations of the present study is the impossibility of comparing the results with the findings of other researchers, as this study has been performed in a circular flume. In addition, this study delivered detailed information on dispersion coefficient on a small scale with low width and there a crucial need to study the advection–dispersion process in high-width channels.

All relevant data are included in the paper or its Supplementary Information.