Turbulence modelling for depth-averaged velocity and boundary shear stress of a dense rigid grass bed open channel

Hydraulic structures in open channels necessitate a thorough understanding of flow characteristics to ensure proper design. Parameters such as mean flow pattern, depth-averaged velocity (DAV) distribution, secondary flow properties, boundary shear stress (BSS), turbulent characteristics and conveyance capacity play a crucial role during the design process. In the study of turbulent flows in open channels, hydraulic engineers often rely on experimental and numerical investigations to understand the various flow components (Guo & Julien 2003; Yang et al. 2004; Sahu et al. 2014; Devi & Khatua (2016); Khuntia et al. (2018a,b); Tajnesaie et al. 2020; Khuntia et al. (2021); Devi et al. (2021); Qasim et al. 2022). Various research has been carried out on the flow and turbulent structure and it influences the riverine habitat (Bornette & Puijalon 2011), water quality (Dosskey et al. 2010), soil erosion and stability (El Bouanani et al. 2022; Ikhsan et al. 2022) and river planform disorderness (Nandi et al. 2022a).

With the advancement of computational techniques, numerical methods have become valuable tools for solving complex hydraulic engineering problems. Even with sophisticated turbulence models, accurately predicting turbulent structures and BSS in free surface flows remains a difficult task (Shnipov 1989; Yang et al. 2004; Yang & Lim 2005; Guo & Julien 2008). Additionally, incorporating vegetation adds complexity to the study of flow-vegetation interaction in open channels. Vegetated flow, which encompasses both physical and biological aspects, has emerged as a significant research area in the field of hydrodynamics (Nepf 1999; Ghani et al. 2019; Nandi et al. 2022b; Wang et al. 2022). In the existing literature, the focus has primarily been on vegetation as roughness elements that result in a decrease in velocity, affecting turbulence parameters (Ree 1958; Chen 1976; Vargas-Luna et al. 2018; Lera et al. 2019). However, it is crucial to recognize that vegetation has a more profound impact on the flow structure due to the generation of coherent structures and redistributes Reynolds shear stress components.

Various approaches have been developed to simulate turbulent flows due to the easiness in numerical modelling. Direct Numerical Simulations (DNS) are comprehensive but computationally intensive models that account for all scales of motion. On the other hand, Reynolds-averaged Navier–Stokes (RANS) models, which assume isotropic turbulence, have become industry standards due to their simplicity. RANS models based on two equation models, such as the k–ε and k–ω models, have been widely used. Reynolds stress models (RSMs) avoid the isotropic assumption and provide a solution for the transport of individual stress terms. Large Eddy Simulation (LES) models simulate most scales by modelling small eddies using sub-grid scale models. Although RSM and LES models tend to yield better results than two equation models, they are typically employed in cases where the application of two equation models is limited, such as in flow separation zones and mixing layers. Previous research by Sofialidis & Prinos (1998) comparing experimental results with a non-linear k–ε turbulence model of low-Reynolds type revealed that the model could not accurately predict the strength of secondary currents observed in experiments. Shiono & Knight (1988) developed the popular analytical model known as the ‘Lateral Distribution Method’, which calculates DAV distributions and later modified it to account for the effects of secondary currents, resulting in the ‘Modified Lateral Distribution Method’. Other researchers, such as Yang et al. (2004), have employed empirical models for secondary flow terms and parabolic distributions for eddy viscosity to simplify RANS equations. However, these methods have encountered difficulties in accurately predicting velocity in rough bed channels (Bonakdari et al. 2008; Guo & Julien 2008; Absi 2011; Lassabatere et al. 2012). Kang & Choi (2006) conducted a comparative study between RSMs and experimental data, demonstrating good agreement between computed and measured results.

With advancements in computational technology, Computational Fluid Dynamics (CFD) techniques have improved the simulation of fluid behaviours in various channel types, rivers and hydraulic structures (Nguyen et al. 2014). Numerical simulations offer cost and time savings while providing detailed information on the flow field with adequate precision (Thanh & Ling-Ling 2015). Sahu et al. (2014) utilized LES and ANFIS models to model velocity distribution and composite friction factor in an asymmetric compound channel, demonstrating the potential for accurate predictions with fewer computational resources using ANFIS. Two equation models, such as the k–ε and k–ω models, have been developed based on the Boussinesq approximation (Qu 2005). The Reynolds stress turbulence model (RSM) has proven effective in capturing three-dimensional flow characteristics and turbulence properties (Asif et al. 2020). Singh & Khatua (2021) and Rezaei & Amiri (2018) predicted the DAV and the BSS distribution in their experimental channels using the standard k–ε model. Tajnesaie et al. (2020) employed four turbulence models to simulate experimental data, comparing secondary flow cells, free-water surface in channels, DAV and shear stress distribution.

However, limited literature exists on the comparative study of various turbulence models for rough open channel flow with grass beds. Therefore, there is a need to analyse the efficiency of different turbulence models to determine their applicability in such conditions. This paper presents an attempt to simulate the experimental results of longitudinal velocity profiles, DAV, BSS and secondary current vectors using seven sub-models from the k–ω and k–ε groups of turbulence models designed for a grass bed open channel. Overall, this paper contributes to the advancement of knowledge in the field of hydraulic engineering by providing valuable insights into the selection and application of turbulence models in rough open channel flow with grass beds.

The structure of this article begins with an introduction section that provides essential background information, literature reviews, establishes the research context and outlines the objectives and novelty of the present study. The subsequent sections describe the experimental set-up, theoretical analysis, modelling approach, analysis of results and conclude with a summary section.

Uniformity of the flow is maintained by installing a tail gate at the downstream end of the flume. To attain the required depth of flow in a steady uniform case, the tailgate was fixed at a specific height. A flow depth of 11 cm is fixed by adjusting the inlet valve. A SonTek Micro 16-MHz acoustic Doppler velocimeter (ADV) is used to measure three-dimensional velocities. The maximum sampling rate of the instrument was adopted, that is 50 Hz. ADV's remote sampling volume is located 5 cm below or above the acoustic transmitter for minimized flow interference. The down looking probe of the instrument is only used to measure the velocities at lower portion of the channel. It is inaccessible for a down looking probe of ADV to measure the upper layer velocity up to 5 cm from the free surface. So, an up-looking probe of ADV is utilized to prevent this limitation by taking the measurements of the upper part. In the present study, one section at X = 7 m from the inlet is chosen to measure the flow parameters where the full boundary layer and secondary flow developments are supposed to be ensured. The experimental data of the same location are also compared with the numerical results. The lateral distribution of velocities was measured at X = 7 m from the entrance (see Figure 1). The velocities at the flume's boundary were then measured using a Preston tube with an outer diameter of 4.77 mm. Table 1 shows the geometry, roughness and hydraulic parameters of the experiment. A post-processing software called WinADV was used to filter the ADV data by presuming a correlation of 70% and an Sound to Noise Ratio (SNR) level of 15 dB.

When using CFD for 3D numerical modelling, the solution of Navier–Stokes (N-S) equations is required, and these equations are based on the laws of mass and momentum conservation inside a flowing fluid. For an incompressible Newtonian fluid turbulent flow, the Navier–Stokes equation is as follows:

To simulate the flow filed in a grass bed channel, different turbulent models have been utilized in the present study. Out of all the available models in the Fluent software, only the two equation models are applied and discussed. The two equation models discussed in this study operates under assumption of isotropic eddy viscosity, the kind of secondary currents generated by the combined effect of geometry and roughness can be sufficiently reproduced by boundary conditions used in this study.

The limitation of the standard k–ε model is that the turbulent fluctuations are isotropic in nature. As a result, in turbulent fluid flows, the k–ε model provides erroneous predictions (Mohammadi & Pironneau 1994).

As the limitations of the standard k–ε model became acknowledged, changes to the model were developed to enhance its efficiency. The renormalization group (RNG) k–ε model and the realizable k–ε model are two of these enhanced k–ε models (Yakhot & Orszag 1986; Shih et al. 1995). The RNG k–ε model has been refined to accommodate for low-Reynolds number implications and fast strained flows. The realizable k–ε model contrasts from the standard k–ε model in that it includes a different turbulent viscosity formulation. In the present study, low-Reynolds number corrections are applied, and wall function is an approximate approach in which semi-empirical equations are used to estimate the flow characteristics near the boundaries. Applications of wall functions are incorporated by using the ‘standard wall functions’ approach in k–ε models.

It is worth noting that the k–ω model contains a low-Reynolds number allowance for near-wall turbulence, which eliminates the need for near-wall modelling, however, wall functions can be implemented when required. Limitations of the k–ε model such as low predictability in curved streamlines and extreme pressure gradients can be overcome by application of the k–ω model. Another advantage is that the epsilon equation and the omega equation are interchangeable as omega and epsilon are related. These advantages laid the stone for developing a more flexible k–ω shear stress transport (SST) model.

The SST k–ω model is a modified version of the standard k–ω model developed by Menter (2002). It combines the k–ω and k–ε turbulence models such that the k–ω is used for the boundary layer's inner region and swaps to the k–ε in free shear flow. A blending function, F₁, activates the Wilcox model near the boundary and the k–ε model in the free stream. This guarantees that the correct model is applied throughout the flow field:

• The k–ω model is ideal for simulating flow in the viscous sub-layer.

• The k–ε model is suitable for forecasting flow behaviour in areas away from the wall.

These properties make the SST k–ω model more precise and dependable for a broader type of flow. The governing equations and auxiliary equations of the SST k–ω model are presented below:

• Turbulent kinetic energy

  

  (13)
• Specific dissipation rate

  

  (14)
• F₁ (Blending function)

  

  (15)
• Kinematic eddy viscosity

  

  (16)
• F₂ (Second blending function)

  

  (17)
• P_k (Production limiter)

  

  (18)

With the goal of turbulence model consolidation, ANSYS introduced a new turbulence model family called generalized k–ω (GEKO). GEKO is a two-equation model based on the k–ω model formulation, but it is adaptable to a wide range of flow scenarios.

To set the attributes of flow field surfaces and fully describe the flow domain's physical features, boundary conditions are necessary on all of the computational domain's boundaries. The effort of selecting a sophisticated turbulence model may fail if the boundary conditions are not applied properly (Bates et al. 2005).

• (a)

  Inlet

A pressure outlet has been applied gauge pressure at outlet is set to zero and the free surface level was kept at 0.11 m.

• (c)

  Wall

For the channel walls, a no slip boundary condition has been adopted, which specifies that the velocity of all fluid particles near to the boundary will assume the velocity of the wall, i.e., zero.

• (d)

  Rough walls

Rough walls are inputted with two boundary conditions related to motion (no slip) and nature of bed. In the present study, the value of equivalent sand grain roughness (k_s) is considered as 0.015 m, i.e., the height of roughness (Khuntia 2020).

• (e)

  Free surface

Free surface is treated as a free slip wall with zero gradients to turbulent kinetic energy and applied to simulation as a symmetry boundary condition.

Detailed experimental and numerical results are presented and discussed in this section. Results of longitudinal velocity profile plots, DAV plots, BSS plots and secondary current vector plots are presented from the discussed turbulent models and compared with experimental results.

From Figure 6(a), it is observed that both the RNG k–ε model and the SST k–ω model are predicting the accurate longitudinal velocity profiles at the centre of the channel as compared to other turbulence models. Because near the centre of the channel, which is far away from the boundary, is free from the wall effect. So, both the advanced models are found to perform well at the middle of the channel. Figure 6(a)–6(c) shows that the SST k–ω model is forecasting the better longitudinal vertical velocity profiles in all the vertical positions of the channel studied as compared to the other turbulence models. The SST k–ω model is a fusion of the standard k–ω model and the k–ε model. A blending function has been used here to combine both the models as one model. Near the wall, the blending function activates the standard k–ω model (which is well suited for simulating the flow near boundaries) and the k–ε model away from the boundaries. The RNG k–ε model also shows good agreement with experimental results at all vertical positions, but it shows less accuracy near the wall region.

From these plots, it is observed that a pair of secondary current vortices is formed as an outcome of channel cross-section and the grass bed roughness. These two pairs of vortices are acting in opposite directions and transporting low momentum fluid from side walls towards the quarter of the channel, where high momentum fluid is present. This confirms the findings form the experimental results of Figure 9. Similar phenomena also have been observed in Figures 7 and 8. However, magnitudes of high velocity are seen at some distance from free surface. Thus, it shows the formation of secondary currents resulting in velocity dip phenomenon. Near the side walls, the intensity is very low, and the secondary currents move upward due to resistance of the vegetated bed. The BSL k–ω model predicted the maximum magnitude of in the range of 0.000182 m/s and for the standard k–ε model, it is 0.000153 m/s. The difference in output by different models is due to their respective model constants and assumptions in modelling shear stress terms. Progress in understanding turbulence-driven secondary currents is limited in literature, 3D modelling approaches may help in a better understanding of this phenomenon.

DAV distribution in an open channel is useful in understanding the flow characteristics of the channel. Additionally, DAV distribution can provide insights into the flow behaviour, such as the presence of eddies, turbulent spots and other features that affect the flow properties. The present research work rigorously evaluates the effectiveness of the k–ω and k–ε models in predicting DAV and BSS distribution by comprehensively applying all the sub-models associated with these models.

From Figures 11 and 15 and Table 3, it is clearly observed that the CES model overpredicts the flow because of the non-incorporation of secondary currents of grass bed in modelling the DAV. The error percentage in calculating the depth-averaged velocities is very high as compared to the other turbulent models. CES is unable to take the effect of grass bed and wall effect of turbulence though it is a quasi-2D model based on SKM through simplified RANS equations. The CES uses the calibrating constants which need to be modified suiting to this turbulent flow on a grass bed channel. However, the numerical 3D modelling tool is successful in matching the pattern of DAV distribution of experimental results. Among all the discussed models, the RNG k–ε model (MAPE = 2.254%, RMSE = 0.004, MPE = −0.108%, MAE = 0.297 and NSE = 0.909) and the SST k–ω model (MAPE = 2.663%, RMSE = 0.005, MPE = 0.193, MAE = 0.357 and NSE = 0.886) provided more accurate results. CES performed the least with MAPE = 12.901%, RMSE = 0.020, MPE = −13.791%, MAE = 1.896 and NSE = −1.014. The reason may be the accurate incorporation of secondary currents due to grass bed by the RNG k–ε model. The MSE of DAV prediction also confirms the efficiency of models. This may be the reason that the RNG k–ε model is a modified form of the standard k–ε model which accounts for the different scales of motion through changes to the production term. The RNG model takes the effects of smaller scale motion due to wall and dense grass bed. Among the turbulent models, the BSL k–ω model performed poor as this model is designed for turbulent flows over flat plates or pipes and generally holds good for isotropic and local equilibrium cases. For the present channel, both the wall and dense grass make unsuitable for better prediction of flow by this model. More precision in results may be achieved by further refining the mesh or by the usage of non-structural grids. From the k–ε and k–ω groups of turbulence models, Tajnesaie et al. (2020) used only the standard k–ε and SST k–ω models for numerical validation of experimental results of Tominaga & Nezu (1991) (T03 case). In Figure 13, it is observed that at the centre of the channel, the experimental depth-averaged profile is having a sudden depression. This may be due to the velocity dip phenomenon. It is also observed that both the models used by Tajnesaie et al. (2020) underpredict the experimental depth-averaged values which is noticed in the present study as well.

From Figures 12 and 16 and Table 4, it is further observed that the SST k–ω model predicted the best results of BSS with MAPE = 3.981%, RMSE = 0.029, MPE = −3.853%, MAE = 2.631 and NSE = 0.960. The CES software is unable to predict the BSS efficiently when compared to the experimental results providing with MAPE = 16.25%, RMSE = 0.105, MPE = −16.433%, MAE = 8.986 and NSE = 0.470. In the prediction of BSS, both the models (SST k–ω and RNG k–ε) also provided the best results with the minimum error value found from the SST k–ω model. Because this model is a modified version of the standard k–ω which combines both the k–ω and k–ε models such that the k–ω model is used for the boundary layer's inner region and swaps to the k–ε in the inner regions of flow, making the SST k–ω better in predicting BSS results as compared to other discussed turbulence models. From Figure 14, it is observed that the numerical BSS results are underpredicted at most of the channel except at the centre and near the wall. Due to the velocity dip phenomenon at the centre of the channel, the experimental BSS plot is seen to have a sudden depression at the centre. But the numerical models did not provide such type of results for the study of Tajnesaie et al. (2020). But, in the present study, all the applied turbulent models and CES over predicted the experimental results of BSS.

The present study aims at simulating the longitudinal velocity profiles, DAV, BSS and secondary flow patterns in a grass bed open channel by the experiments and by the standard turbulence models (two-equation) using ANSYS Fluent. The following conclusions are drawn:

• The SST k–ω model demonstrates superior accuracy in forecasting longitudinal velocity profiles at all vertical positions within the grass bed channel when compared to the other seven simulated models. This improved performance may be attributed to the SST k–ω model's unique combination of the standard k–ω model and the k–ε model. The RNG k–ε model also shows good agreement with experimental results at all vertical positions, but it shows less accuracy near the wall region.

• It is observed that a pair of secondary current vortices is formed as an outcome of channel cross-section and the grass bed roughness. These two pairs of vortices are acting in opposite directions and transporting low momentum fluid from side walls towards the quarter of the channel, where high momentum fluid is present (i.e., confirming the work of other researchers). Secondary current flow fields by all the models are found to provide similar patterns but of different magnitudes which may be due to the different wall effect and model coefficients utilized by the respective modelling approaches.

• In the prediction of lateral distribution of DAV and BSS, both the SST k–ω and RNG k–ε models provided satisfactory results. The RNG k–ε model predicted the DAV with more accuracy as it can account the smaller scale motion due to the combined effect of wall and dense grass bed precisely as compared to other models. Similarly, the SST k–ω model predicted the BSS results more accurately than other models as it is a modified version of the standard k–ω which combines both the k–ω and k–ε models such that the k–ω model is used for the boundary layer's inner region and swaps to the k–ε in the inner regions of flow. This reason makes the SST k–ω model better in predicting BSS results compared to other discussed turbulence models.

• The statistical error analysis shows that the RNG k–ε model demonstrates the best prediction for the DAV results compared to other models (MAPE = 2.254%, RMSE = 0.004, MPE = −0.108%, MAE = 0.297 and NSE = 0.909). For the case of BSS prediction, the SST k–ω model predicted the best (MAPE = 3.981%, RMSE = 0.029, MPE = −3.853%, MAE = 2.631 and NSE = 0.960). It is further observed that the 3D turbulence models provided better results than the quasi-2D model (CES software) in predicting both DAV and BSS. The performance of CES can be improved by incorporating appropriate values of the calibrating coefficients suitable for this grass bed channel condition.

• As there is limited research available to verify the applicability and usability of the RANS-based two-equation turbulence models of k–ω and k–ε groups for a grass bed channel in predicting various flow parameters, the present article will be very helpful to the future researchers in selecting the proper turbulence model for their research work. Future research can be carried out by considering different geometry and roughness conditions (i.e., mobile bed conditions, rigid and flexible vegetations) for better validation of the two-equation turbulence models of k–ω and k–ε groups.

The research was supported by the Department of Science and Technology (DST), Science and Engineering Research Board (SERB), Govt. of India (File No.: CRG/2021/003150 and SERB Qualified Unique Identification Document: SQUID-1969-KK-5771).

This article does not contain any studies with human participants or animals performed by any of the authors.

Informed consent was obtained from all individual participants included in the study.

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

The authors declare there is no conflict.