Hydraulic engineering applications require a good knowledge of turbulent behaviour in non-prismatic channels. This paper aims to predict turbulent behaviour using large-eddy simulation (LES). The model channel has a warped transition. We perform two-phase LES of free-surface flow and validate the results using experimental data and benchmark solution. We discuss rigorous strategies for model set-up, parameter selection and parametric value assignment, including parameters in spectrum synthesiser (SS) and vortex method (VM) for inlet turbulence. The predicted flow displays complex structures due to eddy motions translated from upstream and locally generated by asymmetrical separation in the transition. The history of the flow dynamics may affect the flow development. The predicted velocity, energy spectrum, root-mean-square error, hit-rate and factor-of-two compare well with measurements and benchmark solution. Mapping mean-velocity distribution from experimental data, combined with SS, gives satisfactory inlet condition; alternatively, a 1/7th power-law for the mean-velocity, combined with VM, is acceptable. This paper uses the Okubo–Weiss parameter to delineate instantaneous coherent structures. The LES methods are reliable, efficient and cost-effective. As compared to the simulation of prismatic channels, the flow dynamics in non-prismatic channels exhibit flow separation and turbulence interactions, which increase the flow-complexity, while offering results with crucial practical applications.

  • Key features of turbulence are revealed by two-phase LES of free surface flow in expansion.

  • Effective and practical methods for open-boundary treatments, flow initialisations and turbulence closures are explored.

  • Instantaneous eddy structures are delineated using the Okubo–Weiss parameter.

  • An upstream channel section eight times the flow depth is recommended for LES of flow in expansion.

Open channels need a transition for practical hydraulic engineering applications such as hydro-power generation, water supply, drainage and irrigation. As an example of water supply, water from a lake emerges from a relatively small rectangular concrete channel section and flows through a transition connected to a large trapezoidal earthen channel section. In Figure 1, a warped transition provides the needed connection; the channel is non-prismatic or the channel has an inlet and outlet of different shapes and sizes, in contrast to a prismatic channel (Sahu et al. 2014). To reduce construction costs, real transitions have a short length L and a large divergence angle . This causes flow separation from sidewalls, turbulent disturbances, free-surface fluctuations, poor hydraulic efficiency and high risks of earthen channel erosion.
Figure 1

LES model channel. (a) 3D view; (b) elevation view and (c) top view. The channel dimensions are: ; ; ; and . The green squares mark nine cells: three at , three at 0.47 and three at 0.06. Please refer to the online version of this paper to see this figure in colour: https://dx.doi.org/10.2166/hydro.2023.018.

Figure 1

LES model channel. (a) 3D view; (b) elevation view and (c) top view. The channel dimensions are: ; ; ; and . The green squares mark nine cells: three at , three at 0.47 and three at 0.06. Please refer to the online version of this paper to see this figure in colour: https://dx.doi.org/10.2166/hydro.2023.018.

Close modal

Traditional studies of channel transitions take the empirical approach, which consider only 1D time-averaged free-surface flow (Chow 1959). They ignore 3D instantaneous turbulent features and fail to achieve optimal hydraulic efficiency and minimise erosion risks. The purpose of this study is to predict the turbulent features by means of two-phase large-eddy simulation (LES). The model domain is illustrated in Figure 1. The scope of work includes rigorous validations of LES predictions using newly available experimental data and benchmark solutions. These represent unique and new contributions. In the broad sense, this study demonstrates the two-phase LES approach as a modern approach to future hydraulic design, offering a good complement to the experimental approach.

Unlike the Reynolds-averaged equation models that simulate time-averaged flow field, LES provides instantaneous flow field. Knowledge about the instantaneous flow velocity and the associated peak strength is most relevant to addressing channel erosion risks. Previously, some researchers (e.g., Nashta & Garde 1988) obtained semi-analytical functions for describing transition profiles. Others (Murty Bhallamudi & Hanif Chaudhry 1992; Younus & Hanif Chaudhry 1994; Rahman & Hanif Chaudhry 1997) numerically simulated two-dimensional (2D) depth-averaged flow in expansion. These studies do not produce detailed structures of the flow field or separation features. This is one of the significant gaps addressed by this study.

Compared to other types of transitions (Chow 1959), the warped transition is a superior choice (USCAE 1991). Its sidewalls are flush with both the upstream and downstream channel sections (Li 2022). The seamless joints help reduce turbulent disturbances. A few recent studies (Choi et al. 2022; Zeng & Li 2022) compute instantaneous flows in a straight-wall transition rather than a warped transition.

We address the question of a realistic flow condition at the inlet (Figure 1), with required fluctuations in LES. Periodic conditions often used in LES studies (Xie et al. 2013; Nikora et al. 2019; Ramos et al. 2019) are invalid for a non-prismatic channel. However, inlet conditions may be generated using another method such as precursor method, recycling and rescaling method, or a synthetic method. The 1st method is reportedly the most accurate (Tabor & Baba-Ahmadi 2010), and the 2nd method is effective in producing a fully developed turbulent boundary layer, but both incur high computing costs and can possibly produce a periodicity (Dhamankar et al. 2018). The synthetic methods are a trade-off between computational accuracy and cost. The synthetic methods including a spectral synthesiser (SS) and a vortex method (VM) have been applied in recent LES studies of one-phase air flow (Antoniou et al. 2017).

We also address the question of the possible reflection of surface waves from the outlet (Figure 1) and the influence of initial conditions in LES on the flow development. The results from this study include 3D velocity, vortical motions and first- and second-order statistics of turbulence.

Few applications of SS and VM have been reported in the literature, and little is known about suitable values for parameters as input to the SS and VM. Previously, Liu et al. (2021) applied an SS in their study of flow in a rectangular vegetated channel, and Zeng & Li (2022) applied VM in their investigation of flow choking in a circular partially full pipe. The rectangular channel and circular pipe are prismatic channels. Non-prismatic channels are a new case, where more complicated mechanisms of flow separation, eddy formation and turbulence interaction are involved. This is one important difference of this study from the previous studies. Another important difference of this study is a systematic assessment of proper parameters in the SS and VM and their values. The contributions from this study include improved strategies for computing open channel turbulence and for visualising and analysing the computational results, with relevance to practical applications. The contributions will be discussed later in this paper.

Model equations

We considered an incompressible flow of an air–water mixture, with air above the free surface and water below. The position of the free-surface was tracked using the volume of fluid method. The density, , of the mixture was calculated based on the volume fraction, , as , where the subscripts 1 and 2 refer to air and water, respectively. For a given cell, The viscosity, , of the mixture was calculated as .

LES used filtered equations of momentum and continuity:
formula
(1)
formula
(2)
where is the resolvable-scale filtered velocity component in the -direction ( = 1, 2, 3); t is time; p is the resolvable-scale pressure; is the SGS stress and is the gravitational acceleration. LES computes the exact motion of large eddies and models SGS eddies smaller than the mesh size .

Turbulence closure models

The term accounts for the effect of unresolved velocity fluctuations on the resolved motion. It is related to the resolved strain rate, , through eddy viscosity, , as
formula
(3)
where is the isotropic part of the SGS stress, not modelled but added to the filtered static pressure term. We assess three SGS stress models for : dynamic Smagorinsky-Lilly (DSM) model (Germano et al. 1991); wall-adapting local eddy-viscosity (WALE) model (Nicoud & Ducros 1999) and dynamic kinetic energy (DKE) SGS model (Kim et al. 1997). The 1st and 2nd models predict the behaviour of correctly near a solid surface, without using empirical damping functions, whereas the 3rd model takes the specified filter width, , as the relevant length scale. The standard Smagorinsky SGS stress model has been previously used in LES studies of free surface flow (Cataño-Lopera et al. 2017), where the eddy viscosity coefficient is set as a constant. This is unrealistic for strongly separated flow. Therefore, we do not consider the standard model.

Boundary conditions

We assess three ways to specify the inlet condition for the distribution of mean-inflow velocity :

  • (1)

    1/7th power-law (Table 1, labelled as ‘1/7’);

  • (2)

    1/14th power-law (labelled as ‘1/14’);

  • (3)

    mapped 3D velocity from the laboratory measurements (labelled as ‘EXP’).

Table 1

Boundary condition, initial condition (IC) and hydraulic condition for LES runs

RunITFICSGSNodesNBTReFr
ID(m/s)(L/s)(m)    
T1 0.487 22.6 2.67 1/7th VM ‘L’ DSM 5.13 0.2813 3763707 – 111230 0.33 – – 
T2 VM 7.7 – 111230 – – 
T3 VM 11.6  111230 – – 
T4 SS 5.13  111230 – – 
T5 SS 7.7  111230 – – 
T6 SS 11.6  111230 – – 
T7 0.487 22.6 2.67 1/7th SS ‘L’ DSM 7.7 10a 3763707  111230 0.33 – – 
D1 0.511 23.7b 2.67 1/7th VM ‘L’ DSM 7.7 0.2813 3763707 – 116712 0.34 0.74 0.76 
D2 Yes 116712 0.34 0.72 0.70 
D3 0.602 28.0 8c 1/7th VM ‘L’ DSM 7.7 0.2813 4021768 – 137497 0.40 0.73 0.77 
D4 0.487 22.6 1/14th – 111231 0.33 0.54 0.57 
D5 0.487 22.6 2.67 1/7th SS ‘N’ DSM 7.7 0.2813 3763707 – 111231 0.33 – – 
D6 ‘L’ – 111231 0.33 0.69 0.68 
E1 0.481 22.3 2.67 EXP VM ‘S’ DSM 7.7 0.2813 3763707 – 109860 0.32 0.78 0.77 
E2 SS WALE 2.5 0.2284 – 109860  0.85 0.82 
E3 VM WALE 2.5 0.2284 – 109860  0.82 0.81 
E4 VM WALE 2.5 0.2284 Yes 109860  0.80 0.78 
E5 VM WALE 7.7 0.2813 – 109860  0.77 0.77 
E6 SS WALE 7.7 0.2813 – 109860  0.84 0.81 
E7 VM DKE 7.7 0.2813 – 109860  0.78 0.78 
E8 SS DKE 7.7 0.2813 – 109860  0.87 0.81 
E9 0.481 22.3 8c EXP VM ‘S’ DSM 7.7 0.2813 4021768 – 109860 0.32 0.85 0.84 
P1 0.481 22.3 2.67 INSd – ‘S’ DSM – – 3763707 – 109860 0.32 0.85 0.84 
P2 INSe – – – – 109860 0.32 0.88 0.85 
RunITFICSGSNodesNBTReFr
ID(m/s)(L/s)(m)    
T1 0.487 22.6 2.67 1/7th VM ‘L’ DSM 5.13 0.2813 3763707 – 111230 0.33 – – 
T2 VM 7.7 – 111230 – – 
T3 VM 11.6  111230 – – 
T4 SS 5.13  111230 – – 
T5 SS 7.7  111230 – – 
T6 SS 11.6  111230 – – 
T7 0.487 22.6 2.67 1/7th SS ‘L’ DSM 7.7 10a 3763707  111230 0.33 – – 
D1 0.511 23.7b 2.67 1/7th VM ‘L’ DSM 7.7 0.2813 3763707 – 116712 0.34 0.74 0.76 
D2 Yes 116712 0.34 0.72 0.70 
D3 0.602 28.0 8c 1/7th VM ‘L’ DSM 7.7 0.2813 4021768 – 137497 0.40 0.73 0.77 
D4 0.487 22.6 1/14th – 111231 0.33 0.54 0.57 
D5 0.487 22.6 2.67 1/7th SS ‘N’ DSM 7.7 0.2813 3763707 – 111231 0.33 – – 
D6 ‘L’ – 111231 0.33 0.69 0.68 
E1 0.481 22.3 2.67 EXP VM ‘S’ DSM 7.7 0.2813 3763707 – 109860 0.32 0.78 0.77 
E2 SS WALE 2.5 0.2284 – 109860  0.85 0.82 
E3 VM WALE 2.5 0.2284 – 109860  0.82 0.81 
E4 VM WALE 2.5 0.2284 Yes 109860  0.80 0.78 
E5 VM WALE 7.7 0.2813 – 109860  0.77 0.77 
E6 SS WALE 7.7 0.2813 – 109860  0.84 0.81 
E7 VM DKE 7.7 0.2813 – 109860  0.78 0.78 
E8 SS DKE 7.7 0.2813 – 109860  0.87 0.81 
E9 0.481 22.3 8c EXP VM ‘S’ DSM 7.7 0.2813 4021768 – 109860 0.32 0.85 0.84 
P1 0.481 22.3 2.67 INSd – ‘S’ DSM – – 3763707 – 109860 0.32 0.85 0.84 
P2 INSe – – – – 109860 0.32 0.88 0.85 

All runs have a large Reynolds number (Re = ) and Froude number (Fr) below unity.

aViscosity ratio () is used instead of l or D.

bA larger value is used to compensate lower velocities in sidewall boundary layers.

cThe inlet in Figure 1 is moved to .

dInlet instantaneous velocities are constructed from the LES results for E1.

eInlet instantaneous velocities are constructed from experimental data.

The measurements were made from 70 points in the plane (Figure 1). Inlet turbulent fluctuations (ITFs) are usually needed for the superimposition to . We assess two ways to meet this need: using SS to generate ITF from a total of 100 Fourier harmonics; using 2D VM to superimpose a fluctuating vorticity field from a total of 190 vortices. In both, either the turbulence intensity and the eddy-to-molecular viscosity ratio , and the length scale l, and the hydraulic diameter D, or TKE k and turbulent dissipation rate were supplied. Estimation from the measurements (Ramamurthy et al. 2017) from two planes ( and 0.885) gave ; ; . A Reynolds-averaged simulation suggested .

We assess two ways to specify outlet condition: (1) pressure outlet , with the flow depth equal to 0.2371 m, and (2) pressure outlet, in combination with numerical beach technique (NBT) for damping waves (Park et al. 1999). Other boundary conditions were: on the no-slip bottom and at the no-slip sidewalls. The pressure is zero at the free surface.

LES mesh

The model domain (Figure 1(a)) was discretised into a mesh of hexahedrons. The regions near the bottom, sidewalls were refined using 30 inflation layers. The first layer off a wall had a wall distance , where y is the distance to the wall, and is the friction velocity. The mesh resolves the viscous sublayer. The between-layer spacing varied smoothly. The mesh sizes satisfy Zang's (1991) guidelines and/or Piomelli's (1993) criteria. The growth rate of sizes was kept within 1.2. The mesh had a low skewness; the average value is below 0.17 and the maximum is below 0.65. With low skewness (Chung 2002), the mesh quality is good for maintaining numerical accuracy and stability.

Initial conditions

We assess three different ways to specify the initial conditions:

  • (1)

    The left flank had ; the middle region and right flank had water velocity and , where was uniform at = constant and satisfied (Figure 1(a)). The free surface rose linearly along the transition length, from at to at , where is the flow depth at the outlet; for ; for . was uniform in the -direction. The pressure p was hydrostatic at ; and at . This set of initial conditions is labelled as ‘L’ in Table 1.

  • (2)

    For the entire channel, was uniform at a given cross-section ( = constant) and satisfied . Other remarks are similar to those given in initial condition ‘L’. This set of initial conditions is labelled as ‘N’.

  • (3)

    The ensemble averages of from run D1 (Table 1) were used as initial conditions for runs E1 − E8. This set of initial conditions is labelled as ‘S’.

LES runs

The time step was = 0.01 s for all the LES runs (Table 1), which satisfies the Courant–Friedrichs–Lewy criterion for numerical stability. The model equations were solved in a pressure–velocity coupling manner, using the SIMPLE algorithm. Iterations continued until the convergence criterion of 10−6 or a maximum of 120 iterations per time step was reached. The LES runs were carried out using Ansys FLUENT software (Ansys 2018) and the outputs were visualised and processed using MATLAB (MathWorks 2022) scripts.

Validation of inlet turbulent fluctuations using experimental data

Runs T1 − T3 differed only in turbulence intensity as an input to VM (Table 1). The middle value of was an estimate from experimental data (Ramamurthy et al. 2017). The other input to VM was the hydraulic diameter D. The mean-inflow velocity distribution matched the experimental discharge. Runs T4 − T6 matched the conditions of runs T1 − T3, respectively, but used SS rather than VM. Run T7 matched run T5 but used the viscosity ratio () rather than D.

Commencing from the initial condition ‘L’, runs T1 − T7 allowed a spin-up time period or 26ta, where is the advection time scale , and continued for 2. From output of instantaneous velocity (Equation (1)) over this 2 period, we calculated the ensemble average and the root-mean-square deviation (RMSD) (defined as ), using ensemble size = 250 at a given cell (). Vertical profiles of and for selected () locations for runs T2 and T6 are plotted in Figure 2. The selected locations are well upstream of the transition, free of influence by intense local instability and thus suitable for testing the suitability of ITF. For run T2, the computed profiles correlate well with the experimental data points (Figure 2(a) − 2(e)). The coefficient of correlation was as high as 0.96 (Figure 2(b)) and as low as 0.41 (Figure 2(e)), and the overall average root-mean-square error (RMSE) was as low as for the profile in Figure 2(c) and as high as for the profile in Figure 2(e). The computed profiles plot through the data points (Figure 2(k)–2(o)). Among runs T1–T7 (Table 1), the results for run T2 compare the best with the experimental data.
Figure 2

Vertical profiles of (a) − (j) computed ; (k) − (t) computed at five selected locations marked by triangles in Figure 1(c). Panels (a)–(e) and (k)–(o) are for run T2. Panels (f)–(j) and (p)–(t) are for run T6. Experimental data points (Ramamurthy et al. 2017) are shown for comparison. Please refer to the online version of this paper to see this figure in colour: https://dx.doi.org/10.2166/hydro.2023.018.

Figure 2

Vertical profiles of (a) − (j) computed ; (k) − (t) computed at five selected locations marked by triangles in Figure 1(c). Panels (a)–(e) and (k)–(o) are for run T2. Panels (f)–(j) and (p)–(t) are for run T6. Experimental data points (Ramamurthy et al. 2017) are shown for comparison. Please refer to the online version of this paper to see this figure in colour: https://dx.doi.org/10.2166/hydro.2023.018.

Close modal

Using SS, run T6 produced profiles of similar accuracy (Figure 2(f)–2(j)) as run T2. However, the use of SS surprisingly underpredicted (Figure 2(p)–2(t)). The turbulence intensity was actually large for run T6 than run T2 (Table 1). Among runs T1–T7, run T6 gave the largest underpredictions of . Run T6 has an average for the five selected () profiles. Among the same seven runs, run T3 gave the lowest value (0.67) for the coefficient of correlation, among the five selected profiles. The results for run T7 had similar quality as runs T4 − T6. Ai & Mak (2015) reported severe underpredictions by SS compared to VM in their LES study of air flow over a building, without offering validation. We validate the advantage of VM with turbulence intensity () and hydraulic diameter (), and use them for most of the subsequent runs D1–E9 (Table 1).

Instantaneous velocity

Runs D1–P2 output instantaneous , p and values, following a spin-up period of 56 (long enough to ensure fully developed turbulence). As an example, instantaneous values from three selected cells for run D1 are plotted (Figure 3(a)–3(c)) as a time series. The time was normalised using , which is a time scale related to the fluctuation frequencies of energy-bearing eddies larger than those in the inertial subrange. The value of was determined from the cut-off frequency indicated by the dashed line in Figure 4(a). The inertial subrange is to the right of this dashed line. The maximum value (about 1 s) of was used for normalisation.
Figure 3

(a)–(c) Instantaneous values for run D1 at three selected locations, marked as green squares in Figure 1(c); (d)–(f) cumulative average corresponding to (a)–(c), respectively; (g)–(i) autocorrelation function for time lag k corresponding to (a)–(c), respectively. The coordinates of the locations for panels (a)–(c) are ) = (0.18, 2.44, 0.83), (0.53, 2.51, 0.83), (0.96, 2.66, 0.83). Please refer to the online version of this paper to see this figure in colour: https://dx.doi.org/10.2166/hydro.2023.018.

Figure 3

(a)–(c) Instantaneous values for run D1 at three selected locations, marked as green squares in Figure 1(c); (d)–(f) cumulative average corresponding to (a)–(c), respectively; (g)–(i) autocorrelation function for time lag k corresponding to (a)–(c), respectively. The coordinates of the locations for panels (a)–(c) are ) = (0.18, 2.44, 0.83), (0.53, 2.51, 0.83), (0.96, 2.66, 0.83). Please refer to the online version of this paper to see this figure in colour: https://dx.doi.org/10.2166/hydro.2023.018.

Close modal
Figure 4

Velocity spectral densities of for run D1 (Table 1): (a)–(i) for the nine locations marked as green squares in Figure 1(c). Panels (g)–(i) correspond to Figures 3(a)–3(c), respectively. Please refer to the online version of this paper to see this figure in colour: https://dx.doi.org/10.2166/hydro.2023.018.

Figure 4

Velocity spectral densities of for run D1 (Table 1): (a)–(i) for the nine locations marked as green squares in Figure 1(c). Panels (g)–(i) correspond to Figures 3(a)–3(c), respectively. Please refer to the online version of this paper to see this figure in colour: https://dx.doi.org/10.2166/hydro.2023.018.

Close modal

In Figure 3(a)–3(c), it is shown that the primary flow was stronger than the secondary flow, i.e., and , as expected. The data series of both and had a nearly zero mean. The primary flow had a reverse velocity as large as 0.2. The data series of , and all exhibited significant fluctuations; had larger RMSDs, with a maximum of 0.25.

We examined a near-bed plane and noticed stronger fluctuations near the left sidewall, compared to those at the middle and those near the right sidewall (Figure 1(c), green squares). Near the water surface, near the left sidewall had a mean of , much smaller than at the middle and near the right sidewall, although the left sidewall location is in a narrower section closer to the transition entrance.

Validation of computed turbulence using a benchmark solution

The Kolmogorov −5/3 law is a well-established spectral distribution of turbulent energy in the inertial subrange (Kolmogorov 1941). This study used the law as a benchmark solution for validating LES predictions. We extracted data series of instantaneous velocities over a period of 10 s from the predictions for a series of cells and carried out a Fourier analysis of the data series. Examples of output velocity spectral density, E, from the analysis are plotted in Figure 4. The distributions show a typical energy spectrum for a turbulent flow, the existence of energy containing eddies and a well-fitted −5/3 slope (dashed lines) in the inertial subrange. We conclude that the LES runs properly resolved a major part of turbulence energy and captured the cascading of energy to smaller eddies.

Data samples from the LES output

Averaging instantaneous LES output is needed for practical use of LES results. What is the required minimum length of data series or minimum sample size? We answer this question by evaluating the accumulative average and the autocorrelation function over a time period T elapsed following a spin-up period . Examples of are plotted in Figure 3(d)–3(f). Take as an example = 80 in Figure 3(d). was the average of 1,000 instantaneous values from = 70 to 80. At ≥ 80, changed very little. This is true for other locations in the model channel. Thus, 10 or 10 s is suitable for ensemble averaging, corresponding to sample size = 1,000.

We confirmed the suitability by further calculating . Examples of values are plotted in Figure 3(g)–3(i). If two quantities separated by k are essentially uncorrelated, the selected number of ensembles N is sufficient and suitable for averaging (Box et al. 2015). measures the correlation between the univariate data series at time t and at time , where k increases from 0, … , 10 s. The wavy curves (Figure 3(g)–3(i)) fluctuated in wavelength and amplitude, depending on the direction and location. However, they had some common features at due to a perfect self-correlation; rk decayed with increasing from zero and became bounded by for s. The curves oscillated within a range of small values for . The oscillations might be associated with the Hurst phenomenon (Nordin et al. 1972; Meneveau & Sreenivasan 1991; Dimitriadis et al. 2021). These low values of towards mean that the instantaneous velocities separated by a 10 s time span were uncorrelated. The preceding discussion served the purpose of confirming the suitability of using 10 s for averaging.

Hereinafter, data sampling from the LES output started after the spin-up period (or ) at a sampling frequency and continued until . In summary, the sampling time duration is . This choice is justified on the following grounds: by an order of magnitude; is 8 times the advection time scale ; the cumulative average of velocity supports the choice s; the sufficiency of s is manifested in the autocorrelation function. The ensemble-averaged quantities presented hereinafter use .

Note that in a large water body, turbulent motions can contain slow oscillations at large scales. To capture such oscillations, the sampling time duration needs to be longer. For example, Soulsby (1980) and Walter et al. (2011) used as long as 8–12 min for oceanic motions in order to capture the desired turbulent length scales and obtain quasi-stationary statistics.

Ensemble-averaged flow field

We investigated the distributions of ensemble-averaged horizontal velocity vector, , by selecting three horizontal planes (near-bottom, mid-depth and near-surface). In Figure 5(a), an example of distributions, together with contours of normalised , is shown. For a given location, the sign function Sgn() assigns a plus (minus) sign to the magnitude if (). A positive (negative) magnitude means flow in the positive (negative) -direction. Before the entrance , was zero at the sidewalls and increased rapidly within a small wall-normal distance. After the entrance , the flow separated, mainly along the left sidewall, having one recirculation zone and persisting a long distance. In the transition, the primary flow had a strong core, deflecting from the centreline to the right (to an observer facing downstream). The strongest velocity magnitude was . The primary flow velocity was asymmetrical about the channel centreline.
Figure 5

Horizontal plane, , showing the distributions of ensemble-averaged: (a) horizontal velocity vector for run E1; (b) TKE for D3 and (c) Okubo–Weiss parameter for D3. Please refer to the online version of this paper to see this figure in colour: https://dx.doi.org/10.2166/hydro.2023.018.

Figure 5

Horizontal plane, , showing the distributions of ensemble-averaged: (a) horizontal velocity vector for run E1; (b) TKE for D3 and (c) Okubo–Weiss parameter for D3. Please refer to the online version of this paper to see this figure in colour: https://dx.doi.org/10.2166/hydro.2023.018.

Close modal

The distributions displayed complex 3D flow structures. The near-bottom flow had small pockets of strong velocities in the transition and nearby regions. There was a hysteresis of the occurrence of separation; the further away from the water surface, the further to the transition entrance. The near-surface flow separated in the transition along both sidewalls and triggered eddies in both flanks. The reverse flow had a maximum magnitude of −0.46.

Extending the upstream section length (from for run D1 to 8 for run D3, Table 1) improves the accuracy of first- and second-order statistics of turbulence, as will be demonstrated later. Here, we compared the distributions of ensemble-averaged TKE in the aforementioned three selected horizontal planes. An example of the distributions is shown in Figure 5(b). In the transition, compared to near-surface and near-bottom, the mid-depth had large values, with a maximum of 0.11 (Figure 5(b)). Large pockets of strong occupied the core of recirculation zones in the left flank. Near the surface, many small pockets of strong flanked a strong velocity core, within which the magnitude of was relatively weak. Near the bottom, there were pockets of weak magnitude in the two flanks due to the proximity of channel boundaries.

We quantitatively located turbulent eddies using the Okubo–Weiss parameter , where is the normal strain component, is the shear strain component and is the relative vorticity. Horizontal areas of , where is a certain threshold, are dominated by eddies, with a core in the local minima of . The areas exhibit large vorticity and strong circulation density of the velocity field. Areas of are dominated by strain. Areas of have background vorticity and strain. We chose , where is the standard deviation of values in the flow area. For more details, refer to Zeng & Li (2022).

As an example, contours of normalised for D3 are shown in Figure 5(c). The normalisation used an eddy length scale m. Eddies larger than can be discerned from the computed flow field. For D3, based on and ; for the aforementioned three horizontal planes.

In the transition, a string of strong alternating eddy and strain pockets originated near the transition entrance and developed along the left- and right side of the core of strong velocity (Figure 5(a)). These eddies and strains were coherent turbulence structures and local-scale flow features. The string on the left side traced a meandering path, whose cross-stream amplitude increased as eddy patterns developed towards downstream. At the transition exit, the cross-stream amplitude of the eddy patterns had length scales comparable to . The patterns' development was accompanied by the weakening of local rotations.

After the transition, local rotations persisted longer near the surface than near the bottom and diminished at about . The string on the right side virtually hugged the sidewall, with little cross-stream shifting of eddy positions at the mid-depth (Figure 5(c)). This was not the case near the bottom and near the surface, where the cross-stream amplitude increased slightly.

Validations of ensemble-averaged flow using experimental data

Vertical profiles of ensemble average from nine locations are plotted in Figure 6. The velocities were strong near the surface and decayed towards the bottom. The profiles compared well with experimental data. The computed values correlated well with the experimental values, with the correlation coefficient and a small RMSE. The relative RMSE was as small as 0.01 before the transition, about 0.06 in the transition and below 0.08 after the transition. We accurately predicted strong reverse flow (Figure 6(d)), with a small relative RMSE (equal to 0.07). For a large number of locations (Figure 1(c), the solid dot symbols), the overall average of the RMSE is 0.17. Some of the profiles (Figure 6) show a maximum velocity below the water surface and relatively large variability. The below-surface maximum velocity can be explained by the velocity dip phenomenon (Tominaga et al. 1989). The large variability is due to the transition impact.
Figure 6

Vertical profiles (—) of for E1 at (a)–(i) nine locations, marked by the red multiplication signs in Figure 1(c). The data (+) from Ramamurthy et al. (2017) are shown for validation.

Figure 6

Vertical profiles (—) of for E1 at (a)–(i) nine locations, marked by the red multiplication signs in Figure 1(c). The data (+) from Ramamurthy et al. (2017) are shown for validation.

Close modal
Vertical profiles of from nine locations are shown in Figure 7. The computed values of had an acceptable comparison with experimental values. The relative RMSE was about 0.08 for the profiles before the transition (Figure 7(a)-7(c)). For some locations after the transition, the computed values had relatively large RMSE. For example, in Figure 7(i), relative RMSE was 0.20, but the computed profile showed the same magnitude of fluctuations as the experimental data points. The fluctuations were more profound after the transition (Figure 7(g)7(i)) than before the transition (Figure 7(a)7(c)). This means that the transition contributed to the growth of turbulence and gave rise to very large . This is particularly true near the surface in the transition (Figure 7(d)7(f)). Overall, the LES results are acceptable in terms of second-order turbulence statistics.
Figure 7

Vertical profiles (—) of computed for run E1 at (a)–(i) nine locations, marked by the blue circles in Figure 1(c). The data (+) from Ramamurthy et al. (2017) are shown for validation.

Figure 7

Vertical profiles (—) of computed for run E1 at (a)–(i) nine locations, marked by the blue circles in Figure 1(c). The data (+) from Ramamurthy et al. (2017) are shown for validation.

Close modal
The ensemble-averaged free-surface positions at 30 locations along the channel centreline for runs E1 and D1 are compared with experimental values in Figure 8(a). The relative errors of the computed values for the two runs are below 1.5 and 1.7%. The maximum error occurred near the transition entrance at for run E1 and at the entrance for run D1. Run D1 used the 1/7th law to give inlet mean-flow condition, whereas run E1 used experimental data (Table 1). The computed values correlated well with the experimental values, with the coefficient of correlation . All the runs (Table 1) gave acceptable values of .
Figure 8

(a) Comparison of computed ensemble-averaged free-surface elevations with experimental values at a series of locations. (b) Comparison of computed instantaneous free-surface profiles along the channel centreline at between runs E3 and E4, which are without and with NBT, respectively (Table 1). Please refer to the online version of this paper to see this figure in colour: https://dx.doi.org/10.2166/hydro.2023.018.

Figure 8

(a) Comparison of computed ensemble-averaged free-surface elevations with experimental values at a series of locations. (b) Comparison of computed instantaneous free-surface profiles along the channel centreline at between runs E3 and E4, which are without and with NBT, respectively (Table 1). Please refer to the online version of this paper to see this figure in colour: https://dx.doi.org/10.2166/hydro.2023.018.

Close modal

Discussion of the quality of primary flow predictions

The hit-rate and a factor of two (Okaze et al. 2021) are useful for assessing the quality of alternative LES set-ups, parameters and assigned parametric values. The better the predicted values of a flow variable compared with measured values, the higher the and scores. The quality is acceptable (i.e., the prediction is in reasonable agreement with the measurement, which is the goal) if and . A short description of the metrics is given below. Let , , …, denote n measured values of from experiments, and , , …, denote the corresponding computed values. The validation metrics are expressed as: , where if or and otherwise, is an allowed relative deviation, and is an absolute deviation; , where = 1 if and otherwise. According to Schatzmann et al. (2010), = 0.25. The threshold is 0.20 for and 0.06 for , as estimated following previous studies. Table 1 lists the scores for the runs, which are based on comparisons with 1,184 measured values of streamwise velocities, lumped together from the 11 dotted cross-sections (CSs) in Figure 1(c). The scores are interpreted as follows:

About initial condition, the uniform distribution (for run D5) is unsatisfactory because both scores are out of acceptable range. A plausible reason is that the turbulence development is history-dependent, as demonstrated in Zeng & Li (2022). Once established, the flow tends to hug the left or right wall, depending on the flow history. In the broad sense, the initial condition used in a LES should reflect the flow history of the experiment in question. Either the initial condition ‘L’ or ‘S’ does in this paper.

For inlet mean-flow, a 1/14th power-law (for run D4) is unsatisfactory because the score is too low. A 1/7th law (for D3) achieves acceptable scores. In fact, this law is a good approximation of turbulent boundary layer profiles at high Re, as is the case in this study. Also, the 1/7th law is easier/more efficient to implement than the well-established logarithmic law.

Extending the upstream section length to eight times the depth (for run E9) from 2.67 times (run E1) improves scores. In particular, the extension improves the prediction accuracy of second-order turbulence statistics. It is worth mentioning that the extension increases computational nodes by 7% and computing time by 37% in this study. For similar LES problems of free-surface turbulence, we recommend an adaption length of about eight times the depth, when the inlet conditions reflect the real flow reasonably well. Some previous researchers suggest an adaptation length of three to four times the building height in air flow problems.

The VM (for run D1) slightly improves scores over the SS (for run D6), when the generated ITF is superimposed on a mean-flow profile (given by e.g., a 1/7th law) that is laterally uniform in the absence of secondary flow and eddy motion. A VM generates stronger ITF than a SS, for the same values of input parameters ( and ). The VM is more reliable in building secondary eddies and associated rotations, representing secondary flows of Prandtl's second kind and causing velocity shear. A similar discussion is given in previous LES studies of air flow turbulence (Thordal et al. 2019).

For the first time, this study explores mapping 3D inlet mean-flow from experimental data. In this new case, the SS seems slightly better than the VM, as indicated by comparing run E2 to run E3, run E6 to run E5 and run E8 to run E7. Plausible explanations are that the mapped velocity already contains secondary flow eddies, and it becomes less critical or even redundant to add turbulent eddies from the VM. The same redundancy does not occur in the SS because it yields little influence on secondary velocity.

All the three SGS stresses models: DSM, WALE and DKE (for runs E1, E5 and E7) enjoy a similar degree of success. This is regardless of whether combined with l or with D is used as parameters, as can be seen from comparing run E3 to run E5. The most likely reason is that the mesh for the runs resolves all the important large eddies. The SGS model for smaller eddies has no critical influence on the LES output.

One way to avoid possible wave reflections from the outlet is to implement NBT (for runs D2 and E4). The scores for runs D2 and E4 are close to those for runs D1 and E3. Thus, it is not critical to use NBT. No waves with an amplitude larger than appear in the domain (Figure 9). Smaller ripples exist before the transition and in the transition. The wavelength is about . The ripples in the transition are due to a rapidly expanding width. The steepness is very small . No ripples fluctuate at a time scale comparable to the time step . The free surface overlaps after between runs E3 and E4 (Figure 8(b)). Thus, there are no significant reflections, true as well for runs D1 and D2.
Figure 9

Deviation, , of the instantaneous free surface () at from the ensemble-averaged free surface over a period of for run D3. Please refer to the online version of this paper to see this figure in colour: https://dx.doi.org/10.2166/hydro.2023.018.

Figure 9

Deviation, , of the instantaneous free surface () at from the ensemble-averaged free surface over a period of for run D3. Please refer to the online version of this paper to see this figure in colour: https://dx.doi.org/10.2166/hydro.2023.018.

Close modal

It is constructive to discuss the and scores for CSs in different zones (Figure 1(c)). For the three CSs in zone 1 before the transition (Figure 1(c)), the scores (>0.95) are much better for runs E1 − E8 than runs D1 − D6, due to the fact that runs E1 − E8 map inlet mean-flow from 3D experimental data, preserving the observed flow behaviour before the transition. The scores are not as high just after the transition to some distance downstream (2.49). The influence of the inlet condition on the local flow behaviour weakens and that of separation strengthens, and it is more challenging to accurately predict the behaviour of stronger separated flow.

After the transition, the scores for zones 2, 4 and 6 are good when the upstream section length is eight times the depth (for run D3). The scores change little by using NBT, which is clear from comparing run E4 to run E3 and run D2 to run D1. It is difficult to achieve high scores for zones 4 and 5, where separation is vigorous. The scores barely pass the acceptable threshold for almost all the runs in Table 1. The score is below the acceptable threshold for zone 6 when using 3D experimental data to give inlet mean-flow and the VM to generate fluctuations (for runs E1, E3 − E5 and E7).

The recommended length of an upstream channel section of eight times the flow depth has practical implications. In the studies of Nashta & Garde (1988) and Shettar & Keshava Murthy (1996), the depths are, respectively, 0.059 − 0.12 m and 0.05–0.06 m. The corresponding recommended lengths are 0.47 − 0.96 m and 0.4 − 0.48 m. In Chow (1959), the flume-to-canal expansion has a depth of 1.15 m, the recommended length being 9.2 m. For a concrete siphon expansion with a design depth of 1.25 m (Hinds 1928), the recommended length is 10 m. This study provides some guidelines for the practical choice of channel lengths.

Limitations and strategies for improvement

For those runs where an SS was used to generate inlet velocity fluctuations, the LES methods give underestimations of RMSD in the upstream channel expansion before the expansion, when compared to measured values on a point-to-point basis. Improvement can be achieved by mapping inlet mean-flow from experimental data of 3D velocity, combined with a VM (for runs E1, E5 and E7). The issue is that such data is often not available in many cases. An alternative is to use a 1/7th power-law together with the extended approach section length (for run D3) or a precursor simulation. One limitation is excessively high computing costs. Future work should investigate the influence of varying downstream channel length on LES results.

For improvement, the conditions for run E9 duplicate those for run E1, except that the upstream section length is 3 and the flow field from run D3 provides the initial condition. Run E9 lasts for 64. For the precursor run P1, run E1 continues for 64, with output at each timestep from the plane = 0.89 (Figure 1(c)). The output gives inlet conditions for run P1 over a period of 64. The conditions for run P2 differ from those for run P1 in that the instantaneous inlet velocities are reconstructed as
formula
(4)
where rand is a random function; and are, respectively, the measured mean velocities and standard deviations from the plane . The reconstruction gives 3D velocities at each time over a period of 64. Runs P1 and P2 use instantaneous inlet velocities at each time. Therefore, there is no need to superimpose turbulence. The and scores of RMSD and for runs E9, P1 and P2 are satisfactory for most locations of the LES domain. Note that Equation (4) draws random floating-point numbers from a normal distribution.

New contributions

This study has made unique contributions in two aspects. First, the two-phase flow with a free surface advances from the commonly used rigid-lid approximation. Free-surface variations are important to non-prismatic channel design. Second, the combination of cumulative averaging, autocorrelation and spectrum analysis is an improved method for determining the sample size for ensemble averaging. It avoids trial-and-error and is more efficient and accurate, compared to the traditional method of trying different averaging times.

The Okubo–Weiss parameter is a new alternative to visualise eddy structures. Isosurfaces of represent rotation-dominated regions. The patterns are similar to those using the Q or the criterion. However, this parameter has the advantage that w has a definite threshold . Compared to , w has a physical meaning. For an incompressible flow, if the vertical velocity is asymptotically small. Distributions of 3D instantaneous w at two different timesteps are illustrated in Figure 10, showing strong eddies in the transition. Strong eddy-dominated regions also appear downstream. The parameter clearly delineates instantaneous coherent structures emerging at the transition entrance and evolving downstream. The structures are particularly relevant to channel erosion assessment. The Q-criterion and -criterion involve more than one user-defined threshold value.
Figure 10

Isosurfaces of w for run D3, showing vortex structures at time: (a) 72 and (b) 80. Please refer to the online version of this paper to see this figure in colour: https://dx.doi.org/10.2166/hydro.2023.018.

Figure 10

Isosurfaces of w for run D3, showing vortex structures at time: (a) 72 and (b) 80. Please refer to the online version of this paper to see this figure in colour: https://dx.doi.org/10.2166/hydro.2023.018.

Close modal

In this paper, two-phase LES of turbulent flow is simulated in a non-prismatic open channel. The significance of this work lies in progressively exploring strategies for LES model set-up, selection of proper parameters and parametric value assignments. We validate the LES predictions using experimental data and benchmark solution. The validated strategies can be used to study turbulence under other similar conditions with good accuracy. We show the practical results exerted from this study, which are related to the velocity distributions, spectral distribution of flow energy and hit-rate for validations in other similar studies.

This study reveals 3D complex turbulence structures in the non-prismatic channel, as a result of processes mingled together, involving eddy motions translated from upstream and locally generated by asymmetrical separation. It also reveals the influence of flow history on the turbulence development. One implication is that there may be a need to test a subtle approach to specify the initial conditions in a two-phase LES. This newly recognised challenge has not yet received much attention in the literature.

We show that mapping mean-velocity distribution from measurements, combined with the spectral synthesiser approach to velocity fluctuations, gives satisfactory inlet conditions. An alternative is a 1/7th power-law, combined with the VM. To improve the prediction and control of turbulent behaviour, we recommend adopting a channel length of eight times the depth of the channel.

For the first time, we report instantaneous TKE and coherent turbulence structure in a non-prismatic channel. The practical significance is that the former is of relevance to turbulence control and hydraulic efficiency improvements, whereas the latter helps address channel erosion and instability risks.

This study received financial support from the Natural Sciences and Engineering Research Council of Canada through Discovery Grants held by S.S.L.

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

The authors declare there is no conflict.

Ansys
2018
ANSYS Fluent Theory Guide, Release 19.2
.
ANSYS, Inc
,
Canonsburg, Pennsylvania
.
Antoniou
N.
,
Montazeri
H.
,
Wigo
H.
,
Neophytou
M. K. A.
,
Blocken
B.
&
Sandberg
M.
2017
CFD and wind-tunnel analysis of outdoor ventilation in a real compact heterogeneous urban area: evaluation using ‘air delay’
.
Building and Environment
126
,
355
372
.
Box
G. E.
,
Jenkins
G. M.
,
Reinsel
G. C.
&
Ljung
G. M.
2015
Time Series Analysis: Forecasting and Control
.
John Wiley & Sons
,
Hoboken, New Jersey
.
Cataño-Lopera
Y. A.
,
Landy
B. J.
&
García
M. H.
2017
Unstable flow structure around partially buried objects on a simulated river bed
.
Journal of Hydroinformatics
19
(
1
),
31
46
.
Choi
B. H.
,
Anand
N. K.
,
Hassan
Y. A.
&
Sabharwall
P.
2022
Large eddy simulation of flow through an axisymmetric sudden expansion
.
Physics of Fluids
34
(
6
),
065117
.
Chow
V. T.
1959
Open-Channel Hydraulics
.
McGraw-Hill
,
Caldwell, New Jersey
.
Chung
T. J.
2002
Computational Fluid Dynamics
.
Cambridge University Press
,
Cambridge
,
UK
.
Dhamankar
N. S.
,
Blaisdell
G. A.
&
Lyrintzis
A. S.
2018
Overview of turbulent inflow boundary conditions for large-eddy simulations
.
AIAA Journal
56
(
4
),
1317
1334
.
Germano
M.
,
Piomelli
U.
,
Moin
P.
&
Cabot
W. H.
1991
A dynamic subgrid-scale eddy viscosity model
.
Physics of Fluids A: Fluid Dynamics
3
(
7
),
1760
1765
.
Hinds
J.
1928
The hydraulic design of flume and siphon transitions
.
Transactions of the American Society of Civil Engineers
92
(
1
),
1423
1459
.
Kim
W. W.
,
Menon
S.
,
Kim
W. W.
&
Menon
S.
1997
Application of the localized dynamic subgrid-scale model to turbulent wall-bounded flows
. In:
Meeting Papers on Disc of the 35th Aerospace Sciences Meeting and Exhibit
,
6–9 January 1997
,
Reno, NV, USA.
AIAA Paper 97-0210. https://doi.org/10.2514/6.1997-210
.
Kolmogorov
A. N.
1941
Dissipation of energy in the locally isotropic turbulence
. In:
Doklady Akademii Nauk SSSR
32
,
19
21
.
Li
S. S.
2022
Harmonic function for 3D warped transition geometry and its practical use
.
Journal of Irrigation and Drainage Engineering
148
(
6
),
06022002
.
Liu
X. D.
,
Liu
Z. Q.
,
Tang
L. C.
,
Han
Y.
,
Chen
J.
&
Yang
S. Q.
2021
Analysis of vegetation resistance based on two typical distribution types in ecological channel
.
Ecological Engineering
169
,
106325
.
MathWorks
2022
MATLAB Version: 9.13.0 (R2022b)
.
MathWorks Inc
,
Natick, Massachusetts
.
Meneveau
C.
&
Sreenivasan
K. R.
1991
The multifractal nature of turbulent energy dissipation
.
Journal of Fluid Mechanics
224
,
429
484
.
Murty Bhallamudi
S.
&
Hanif Chaudhry
M.
1992
Computation of flows in open-channel transitions
.
Journal of Hydraulic Research
30
(
1
),
77
93
.
Nashta
C. F.
&
Garde
R. J.
1988
Subcritical flow in rigid-bed open channel expansions
.
Journal of Hydraulic Research
26
(
1
),
49
65
.
Nicoud
F.
&
Ducros
F.
1999
Subgrid-scale stress modelling based on the square of the velocity gradient tensor
.
Flow, Turbulence and Combustion
62
(
3
),
183
200
.
Nikora
V. I.
,
Stoesser
T.
,
Cameron
S. M.
,
Stewart
M.
,
Papadopoulos
K.
,
Ouro
P.
,
McSherry
R.
,
Zampiron
A.
,
Marusic
I.
&
Falconer
R. A.
2019
Friction factor decomposition for rough-wall flows: theoretical background and application to open-channel flows
.
Journal of Fluid Mechanics
872
,
626
664
.
Nordin
C. F.
,
McQuivey
R. S.
&
Mejia
J. M.
1972
Hurst phenomenon in turbulence
.
Water Resources Research
8
(
6
),
1480
1486
.
Okaze
T.
,
Kikumoto
H.
,
Ono
H.
,
Imano
M.
,
Ikegaya
N.
,
Hasama
T.
,
Nakao
K.
,
Kishida
T.
,
Tabata
Y.
,
Nakajima
K.
&
Yoshie
R.
2021
Large-eddy simulation of flow around an isolated building: a step-by-step analysis of influencing factors on turbulent statistics
.
Building and Environment
202
,
108021
.
Park
J. C.
,
Kim
M. H.
&
Miyata
H.
1999
Fully non-linear free-surface simulations by a 3D viscous numerical wave tank
.
International Journal for Numerical Methods in Fluids
29
(
6
),
685
703
.
Piomelli
U.
1993
High Reynolds number calculations using the dynamic subgrid-scale stress model
.
Physics of Fluids A: Fluid Dynamics
5
(
6
),
1484
1490
.
Rahman
M.
&
Hanif Chaudhry
M.
1997
Computation of flow in open-channel transitions
.
Journal of Hydraulic Research
35
(
2
),
243
256
.
Ramamurthy
A. S.
,
Thapa
D. R.
&
Li
S. S.
2017
Experimental study of flow past a warped transition
.
Journal of Irrigation and Drainage Engineering
143
(
8
),
04017022
.
Ramos
P. X.
,
Schindfessel
L.
,
Pêgo
J. P.
&
De Mulder
T.
2019
Flat vs. curved rigid-lid LES computations of an open-channel confluence
.
Journal of Hydroinformatics
21
(
2
),
318
334
.
Sahu
M.
,
Mahapatra
S. S.
,
Biswal
K. C.
&
Khatua
K. K.
2014
Prediction of flow resistance in a compound open channel
.
Journal of Hydroinformatics
16
(
1
),
19
32
.
Schatzmann
M.
,
Olesen
H.
&
Franke
J.
2010
COST 732 Model Evaluation Case Studies: Approach and Results
.
Meteorological Inst.
Hamburg, Germany
.
Shettar
A. S.
&
Keshava Murthy
K.
1996
A numerical study of division of flow in open channels
.
Journal of Hydraulic Research
34
(
5
),
651
675
.
Soulsby
R. L.
1980
Selecting record length and digitization rate for near-bed turbulence measurement
.
Journal of Physical Oceanography
10
(
2
),
208
219
.
Tabor
G. R.
&
Baba-Ahmadi
M. H.
2010
Inlet conditions for large eddy simulation: a review
.
Computers & Fluids
39
(
4
),
553
567
.
Thordal
M. S.
,
Bennetsen
J. C.
&
Koss
H. H. H.
2019
Review for practical application of CFD for the determination of wind load on high-rise buildings
.
Journal of Wind Engineering and Industrial Aerodynamics
186
,
155
168
.
Tominaga
A.
,
Nezu
I.
,
Ezaki
K.
&
Nakagawa
H.
1989
Three-dimensional turbulent structure in straight open channel flows
.
Journal of Hydraulic Research
27
(
1
),
149
173
.
USACE
1991
Hydraulic Design of Flood Control Channels
.
Engineer Manual No. 1110-2-1601
.
Department of the Army Corps of Engineers
,
Washington, DC
,
USA
.
Walter
R. K.
,
Nidzieko
N. J.
&
Monismith
S. G.
2011
Similarity scaling of turbulence spectra and cospectra in a shallow tidal flow
.
Journal of Geophysical Research: Oceans
116
,
C10019
.
Xie
Z.
,
Lin
B.
&
Falconer
R. A.
2013
Large-eddy simulation of the turbulent structure in compound open-channel flows
.
Advances in Water Resources
53
,
66
75
.
Younus
M.
&
Hanif Chaudhry
M.
1994
A depth-averaged turbulence model for the computation of free-surface flow
.
Journal of Hydraulic Research
32
(
3
),
415
444
.
Zang
T. A.
1991
Numerical simulation of the dynamics of turbulent boundary layers: perspectives of a transition simulator
.
Philosophical Transactions of the Royal Society of London. Series A: Physical and Engineering Sciences
336
(
1641
),
95
102
.
This is an Open Access article distributed under the terms of the Creative Commons Attribution Licence (CC BY 4.0), which permits copying, adaptation and redistribution, provided the original work is properly cited (http://creativecommons.org/licenses/by/4.0/).

Supplementary data