Detached Eddy Simulation Investigation of Turbulence at a Circular Pier with Scour Hole

Local scour at a pier is attributable to the forces exerted on the bed by the complex, highly three-dimensional (3D) and unsteady flow field generated by the pier. The flow field is marked by turbulence structures with a wide range of scales, and by a pronounced downflow at the pier’s leading edge. The present paper presents insights into the flow field at a pier with a developed scour hole formed in a sand bed at a relatively high pier Reynolds number, and discusses the linkages between flow field and sediment entrainment and transport from the scour hole. These linkages are relevant for piers with relatively large scour holes corresponding to the later stages of the scour process. The complexity of the flow field and the presence of a large scour hole hamper the use of existing flow instrumentation to obtain detailed insights into the flow and turbulence structure.

To date, most of the focus on the scour flow field has been on the horseshoe vortex (HV) system forming around the upstream base of the pier. It is thought to be the main mechanism responsible for the growth of the scour hole past the initial stages of the scour process (e.g., Dargahi 1990; Melville and Coleman 2000). However, the HV is part of a system of turbulence structures that together with the downflow and the flow acceleration around the flanks of the pier are the erosive flow mechanisms of primary importance. The turbulence structures [e.g., necklace vortices of HV system, eddies shed in the separated shear layers (SSLs), large-scale wake rollers] are not isolated from each other. They intrinsically connect within the flow field. Thus, it is not enough to focus on one turbulence structure, notably the well-known HV system to understand the physics of scour processes around the pier.

To understand pier scour requires quantitative and qualitative knowledge of the flow field formed at a pier during the scouring process and of the dynamics of the large-scale structures that are responsible for entrainment and transport of most sediment particles. This need is reflected by the numerous experimental studies of flow fields around bridge piers conducted especially over the last couple of decades; many of these studies are summarized in Melville and Coleman (2000) and Sumer and Fredsoe (2002). Dargahi (1990) documented the changes in the scour mechanism from flat to equilibrium scour conditions at a circular pier ( ${\mathsf{R}}_D=UD/\nu =4\times {10}^4$ , $U$ =mean incoming flow velocity; $D$ =diameter or equivalent width of the pier; and $\nu$ =molecular kinematic viscosity) placed in a loose bed channel. In particular, he observed no significant change in the structure of the HV system past the initial stages of the formation of the scour hole. This is a significant result for the relevance of the present investigation that focuses on the dynamics of coherent structures for the case of a scour hole at equilibrium scour conditions. 3D velocity distributions inside an equilibrium scour hole at a circular pier were measured by Dey et al. (1995) using a five-hole pitot tube. Measurements of the mean-flow fields and turbulence around a circular pier with equilibrium scour bed at ${\mathsf{R}}_D\sim {10}^5$ were reported by Graf and Istiarto (2002) in several polar planes $\left(0°<\varphi <180°\right)$ . The bed shear stress distribution was inferred from the velocity measurements (extrapolation toward the bed of the estimated Reynolds stresses at stations situated close to the bed surface). The study by Roulund et al. (2005) provides velocity measurements in the symmetry plane upstream and downstream of the pier using a laser Doppler anemometer and estimations of the bed shear stress around the pier using a hot film probe for a circular bridge pier on a flat bed at ${\mathsf{R}}_D\sim {10}^5$ . Dey and Raikar (2007) investigated the turbulent flow fields around circular piers $\left({\mathsf{R}}_D\sim 4\times {10}^4\right)$ within a developing scour hole having depths of 0.25, 0.5, 0.75, and 1.0 times the depth at equilibrium scour conditions based on acoustic Doppler velocimeter measurements. The characteristics of the HV system with the development of the scour hole were discussed based on the observed similarity with the velocity and turbulence characteristic scales. Unger and Hager (2007) performed a particle image velocimetry (PIV) study of the flow around a circular pier over a wide range of Reynolds numbers $\left(5\times {10}^4<{\mathsf{R}}_D<3.5\times {10}^6\right)$ . The interest was in the visualization of the flow patterns in the vertical symmetry plane and in several horizontal planes as the scour process evolves toward equilibrium conditions. The downflow and the HV system were identified as the main reasons for the growth of the scour hole. They also confirmed the fact that the main scour activity is concentrated in the region situated in front of the pier where the coherence of the HV system is highest. Their experimental observations confirmed the variation in the position of the boundary layer separation point on the pier from $\varphi \sim 90°$ away from the bed to $\varphi \sim 180°$ close to the bed, consistent with the large eddy simulation (LES) results of Kirkil et al. (2008). Several of these investigations found that the structure of the HV system varies little past the initial stages of the scour process.

The highly unsteady 3D flow around a pier around is hard to visualize, let alone measure, in their entirety. Presently, it is not possible, by means of laboratory experiment, to visualize the entire instantaneous flow field around a bridge pier or to obtain accurate information on the flow and turbulence statistics very close to the bed (e.g., distributions of the mean and instantaneous bed shear stress and pressure RMS fluctuations). Simulations using eddy-resolving techniques and sufficiently fine meshes to integrate the equations through the viscous sublayer allow us to obtain this information. Numerical models using eddy-resolving techniques like LES and hybrid RANS-LES techniques (e.g., Rodi 1997; Spalart 2000; Tokyay and Constantinescu 2006; McCoy et al. 2008) have been shown to be much more successful to predict massively separated flows and flows dominated by unsteady large-scale coherent structures compared to the classical Reynolds averaged Navier-Stokes (RANS) approach. Also, the knowledge of the mean-flow and turbulence statistics only is not enough to allow investigating the role played by the energetically important coherent structures in the flow. Even if the evolution of the scoured bed was not part of the solution, as the time scales of the large-scale coherent structures are much smaller than the time scales associated with bathymetry evolution, information from these simulations can help to better understand the scouring mechanism at different stages of the erosion and deposition process around the pier. Together with PIV based studies (e.g., Unger and Hager 2007), high-resolution eddy-resolving numerical simulations are the best approach to elucidate fundamental aspects of momentum and sediment transport mechanisms around bridge piers.

Ideally, sufficiently resolved LES simulations without wall functions and employing Smagorinsky subgrid-scale (SGS) models based on the flow physics like the dynamic Smagorinsky model (e.g., Bou-Zeid et al. 2005) should be used to investigate such complex flows. Kirkil et al. (2008) reported such a study for a pier with scour hole at a relatively low Reynolds number $\left({\mathsf{R}}_D=16,000\right)$ . A similar LES investigation was conducted by Koken and Constantinescu (2008a,b) for abutments and by McCoy et al. (2007) for emerged and submerged groynes exposed to the flow. One important question is to what extent the findings from these LES studies and, for that matter, the findings from detailed flume studies conducted at relatively low Reynolds numbers $\left({\mathsf{R}}_D\sim {10}^4\right)$ are relevant for practical applications. At very large Reynolds numbers $\left({\mathsf{R}}_D>{10}^5\right)$ , the mesh and time-step requirements make the use of LES without wall functions not feasible. To investigate Reynolds-number-induced scale effects using eddy-resolving techniques, LES with wall functions or hybrid RANS-LES methods are required. Hybrid zonal (e.g., Wang and Moin 2002) and nonzonal RANS-LES methods can resolve the dynamically most important eddies in the flow and can be used to study flow past bridge piers at ${\mathsf{R}}_D$ closer to those encountered in rivers.

Nonzonal RANS-LES methods like detached eddy simulation (DES) use the same base turbulence model in the RANS and LES regions and no special treatment is required to match the solutions at the boundary between the LES and RANS regions. The model used in the LES region is a RANS model in which the turbulence length scale is modified to allow the energy cascade to the small scales like in LES. For internal flows (e.g., flow past surface-mounted bodies or bottom cavities in a channel) the predictions of DES, including the dynamics of the large-scale coherent structures, are dependent on the presence, or not, of resolved turbulence in the incoming flow (e.g., Chang et al. 2007). The problem is that it is not easy to generate a flow with the proper turbulence content in a channel at a high Reynolds number. Basically all experimental investigations of local scour at bridge piers are conducted with an incoming turbulent channel flow. Thus, it is important that the inflow conditions in the numerical simulations are as close as possible to those present in experiments or in the field.

Kirkil et al. (2008) reviewed recent RANS investigations of the flow past circular piers with fixed or evolving bed. Their study also summarizes prior studies using LES and hybrid RANS-LES methods for related flows (e.g., flow past abutments, piers of complex shape, and other obstacles mounted on flat beds). LES of the flow past a circular pier with flat and deformed bed was performed by Tseng et al. (2000) and Choi and Yang (2002) using relatively coarse meshes (largest simulations used 250,000 cells), wall functions and a constant coefficient SGS model. The incoming flow did not contain any resolved turbulence. Kirkil et al. (2008) were the first to use sufficiently well-resolved LES (mesh contained close to 4 million cells and ${\mathsf{R}}_D=16,000$ ) to study, besides the mean-flow and turbulence statistics, the dynamics of the large-scale coherent structures in the flow past a circular pier with equilibrium scour hole in a relatively shallow channel. The simulations were conducted using a nondissipative energy-conserving code (Mahesh et al. 2004) and the dynamic Smagorinsky model was used to estimate the SGS viscosity. The incoming flow was fully turbulent and contained realistic fluctuations obtained from a preliminary channel LES simulation. The simulation captured the presence of bimodal aperiodic oscillations inside the HV system region, as well as the associated increase in the resolved turbulent kinetic energy $\left(k\right)$ and pressure RMS fluctuations $\left(\overline{p^{\prime 2}}\right)$ . The study showed, for the first time, that bimodal oscillations of the necklace vortices are present in the later stages of the scouring process and, thus, are a general characteristic of bridge pier flows. Dye visualizations and large-scale PIV were used to confirm some of the flow phenomena identified based on analysis of the LES flow fields. Several mechanisms that can explain the growth of the scour hole laterally and behind the pier were identified. DES predictions of flow past wing-shaped cylinders mounted on a flat surface at ${\mathsf{R}}_D=1.15\times {10}^5$ were reported by Paik et al. (2007) using steady inflow conditions (no resolved inflow turbulence). DES predictions of flow past a vertical spur dike (bridge abutment) in a flat-bed channel at ${\mathsf{R}}_D=5\times {10}^5$ were reported Koken and Constantinescu (2008c, 2009) using inflow conditions that contained resolved fluctuations. A similar investigation was conducted by Kirkil and Constantinescu (2009) for a high aspect ratio rectangular cylinder at high angle of attack in a flat-bed channel at ${\mathsf{R}}_D=2.5\times {10}^5$ .

Development of practical design approaches for predicting scour depths at bridges has been hampered by inadequate knowledge about, and formulation of, important component processes at play during scour. A better understanding of the turbulence characteristics, and in particular of the role of the large-scale coherent structures, should lead to scour prediction methods that incorporate more physics and predict more accurately the scour depth for a broader range of flow characteristics (Ettema et al. 2006). Though not always the case [e.g., scour prediction methods proposed by Melville and Coleman (2000) and by Oliveto and Hager (2002a,b)], some scour prediction methods are calibrated using mostly data from relatively low-Reynolds-number flume studies. Thus, information on scale effects associated with an increase in the Reynolds number obtained from either simulations or experiments should be useful toward developing expressions for correction parameters in these methods. At a more fundamental level, there is a need to understand if the turbulence structure, and thus the sediment entrainment and transport mechanisms, change significantly between the Reynolds numbers at which many flume studies are conducted and Reynolds numbers $\left({\mathsf{R}}_D>{10}^5\right)$ above the range that is known to be subject to significant scale effects in other massively separated turbulent flows [e.g., see Williamson (1996), for flow past infinitely long cylinders].

The DES reported here investigates flow past a circular pier with scour hole corresponding to conditions close to equilibrium scour for ${\mathsf{R}}_D=2.06\times {10}^5$ . The scour hole dimensions were determined from a flume experiment. Following a brief description of the numerical method, DES model and simulation setup, the dynamics of the large-scale coherent structures, the mean-flow and turbulence statistics are analyzed. The main goal of the present paper is to discuss the dynamics of these structures for the case of a cylinder with a developed scour hole, representative of the later stages of the scour process. The dynamics of these structures is expected to change significantly between the initial stages when the HV system is much more unstable (e.g., Kirkil and Constantinescu 2007) and the advanced stages when the scour hole stabilizes the HV system and modifies its interactions with the SSLs. On the other hand, the dynamics is expected to be similar once the scour hole is large enough to stabilize the HV system. Scale effects are discussed based on comparison of DES with results from an LES simulation at ${\mathsf{R}}_D=16,000$ (Kirkil et al. 2008). In both simulations the incoming flow contained resolved fluctuations and the ratio between the channel depth away from the pier $\left(H\right)$ and the diameter of the pier $\left(D\right)$ was the same $\left(H/D=1.12\right)$ . One should point out that, as opposed to the flat-bed case, differences in the two solutions are not due solely to the difference in the Reynolds number but also to differences in the geometry of the scour hole, which is deeper in the high Reynolds number simulation. In the ensuing sections, the DES solution at ${\mathsf{R}}_D=2.06\times {10}^5$ is termed the HR solution while the LES solution at ${\mathsf{R}}_D=16,000$ (Kirkil et al. 2008) is termed the LR solution.

The present study used the Spalart-Allmaras (SA) RANS model as the base model in DES (e.g., Spalart 2000). The SA model is based on a transport equation for the modified eddy viscosity, $\tilde{\nu }$ . The transport equation for $\tilde{\nu }$ iswhere $S=\text{magnitude}$ of the vorticity; $u^j=\text{contravariant}$ resolved velocity; $t=\text{time}$ ; $d=\text{distance}$ to the closest wall; and ${\xi }^j=\text{curvilinear}$ coordinate in the $j$ direction. The other variables and parameters areThe eddy (SGS) viscosity ${\nu }_t$ is obtained fromwhereTo account for roughness effects the distance to the (rough) wall, $d$ , is redefined [see also Spalart (2000)] aswhere $d_{\text{min}}=\text{distance}$ to the closest wall and $k_s=\text{equivalent}$ roughness height (e.g., this term can include the form roughness due to presence of ripples or small dunes that are not resolved by the grid). For smooth walls, $k_s$ is taken equal to zero. The model constants in the aforementioned equations are: $C_{b1}=0.135$ , $C_{b2}=0.622$ , $\sigma =0.67$ , $\kappa =0.41$ , $C_{v1}=7.1$ , $C_{w2}=0.3$ , $C_{w3}=2.0$ , and $C_{w1}=C_{b1}/{\kappa }^2+\left(1+C_{b2}\right)/\sigma$ .

The SA version of DES is obtained by replacing the turbulence length scale $d$ with $l_{\text{DES}}$ which is defined as $l_{\text{DES}}=\text{min}\left(d,C_{\text{DES}}\Delta \right)$ where the model parameter $C_{\text{DES}}$ is equal to 0.65 and $\Delta$ is a measure of the local grid size. When the production and destruction terms of the model are balanced, the length scale in the LES regions $l_{\text{DES}}=C_{\text{DES}}\Delta$ becomes proportional to the local grid size and yields an eddy viscosity proportional to the mean rate of strain and ${\Delta }^2$ as in LES with a Smagorinsky model, which allows the energy cascade down to grid size. The eddy viscosity predicted by DES in the LES regions goes to zero if the local grid size decreases to zero as in classical LES.

A general description of the DES code used in the present work is given in Constantinescu and Squires (2004). The 3D incompressible Navier-Stokes equations are integrated using a fully implicit fractional-step method. The governing equations are transformed to generalized curvilinear coordinates on a nonstaggered grid. The discrete momentum (predictor step) and turbulence model equations are integrated in pseudotime using alternate direction implicit approximate factorization scheme. Time integration is done using a double time stepping algorithm and local time stepping is used to accelerate the convergence at each physical time step. Source terms in the turbulence model equations are treated implicitly. The time discretization is second-order accurate. Detailed grid sensitivity and validation studies for the flow past spheres and for flow and contaminant transport in a channel with a bottom cavity are discussed in Constantinescu and Squires (2003,2004), Constantinescu et al. (2002,2003) and in Chang et al. (2007). In particular, the study of Chang et al. (2007) showed that the agreement between highly resolved LES (Chang et al. 2006) and DES on a much coarser mesh improved significantly when the inflow contained turbulence fluctuations obtained from a preliminary channel flow calculation. This is the approach adopted in the present study.

Experiments were conducted in a flume at Iowa Institute of Hydraulic Research (IIHR). Its test section is 20 m long, $W=3\text{m}$ wide and 2 m deep. The diameter of the circular pier was $D=0.46\text{m}$ , the flow depth was $H=0.53\text{m}$ and the mean incoming channel velocity was $U=0.45\text{m}/\text{s}$ . The main nondimensional flow parameters were ${\mathsf{R}}_D=2.06\times {10}^5$ $\left({\mathsf{R}}_H=2.4\times {10}^5\right)$ and ${\mathsf{F}}_H=U/\sqrt{gH}\sim 0.2$ .

A uniform layer of sediment with a median diameter $d_{50}=1.05\text{mm}$ was placed on the bottom of the flume. The thickness of the sand layer was $t_s=1.0\text{m}$ . The critical friction velocity was $u_{\tau c0}/U=0.056$ , where $u_{\tau c0}$ was estimated based on the Shields diagram. The experiment was conducted under clearwater scour conditions and was run until the scour hole growth was negligible. This process took close to 5 days. The discharge was slowly reduced such that the bathymetry changes during the reduction of the flow discharge were small. The measured decrease of the maximum scour depth during the reduction of discharge and dewatering of the flume was less than 1.5%.

The equilibrium bathymetry is shown in Fig. 1(a). The maximum flow depth in the scour hole was $H_S=2.09D$ corresponding to a scour depth of $0.97D$ . The scour depth is significantly larger than that $\left(0.68D\right)$ present in the experiment conducted at ${\mathsf{R}}_D=16,000$ (Kirkil et al. 2008). This larger value of $H_s$ relative to $H$ is closer to the values generally observed at circular bridge piers. Table 1 compares the main hydraulic and geometric parameters in the experiments performed to obtain the bathymetry in simulations LR and HR.

Bathymetry measurements around the pier were taken at 0.1-m intervals, both in the longitudinal and transversal directions. The bathymetry data were then digitized for the purpose of generating the scoured bed surface used in the numerical simulation [Fig. 1(b)]. The measured bathymetry is essentially symmetric with respect to the $\varphi =0°$ (the polar angle $\varphi$ is measured starting at the symmetry plane upstream of the pier) plane. The scour and deposition patterns are similar to those observed to form in previous experimental studies of scour around piers of circular section (e.g., Dargahi 1990; Graf and Istiarto 2002; Kirkil et al. 2008). The slope of the bed inside this part of the scour hole was around 32°. The scour hole is not as deep behind the pier but it extends to larger distances away from the pier. Several small-scale ripples are observed inside the scoured region downstream of the pier. In the central region behind the pier $\left(x/D\sim 5\right)$ , a deposition region in the form of a dune is present. The axis of the deposition dune is streamwise oriented. The deposition dune is situated farther away from the pier compared to the lower ${\mathsf{R}}_D$ experiment of Kirkil et al. (2008).

The computational domain [Fig. 1(b)] extends $7D$ upstream of the pier axis and $20D$ downstream of it. The domain width is $6.67D$ , the same as that of the flume used to perform the experiment. The depth of the channel in the flat-bed region is $H=1.12D$ . Though this ratio seems small compared to the usual ratios of $H/D$ encountered in practice at circular piers, the dynamics of the HV system within the scour hole is not expected to be significantly affected by the flow shallowness if the flow structures that may form in the upstream stagnation region close to the free surface (e.g., surface roller) do not interact directly with the HV system. Our experiments showed that this is the case once $H>0.7D$ . Still, the relatively shallow flow conditions considered in the present investigation may have a strong influence on some aspects of the flow (e.g., near-wake region). Reducing the flow shallowness in our test case would have required an increase of the total number of computational mesh points needed to sufficiently resolve the flow proportional to the increase in the ratio $H/D$ .

The computational domain contains about 9.2 million cells $\left(480\times 240\times 80\right)$ . The mesh in several representative sections is shown in Fig. 2. The larger number of grid points compared to that used in the LES simulation at ${\mathsf{R}}_D=16,000$ (Kirkil et al. 2008) was needed to resolve the attached boundary layers at ${\mathsf{R}}_D=2.06\times {10}^5$ and to capture the dynamically important eddies forming around the pier. The increase in the mesh size in DES is much smaller than the increase required by a sufficiently well-resolved LES simulation at the same ${\mathsf{R}}_D$ because in DES the thin attached boundary layers are resolved using RANS. The grid spacing in the wall-normal direction was close to 0.25 wall units (assuming the nondimensional friction velocity $u_{\tau }/U$ is 0.05) on all the wall surfaces.

The lateral boundaries corresponding to the position of the flume sidewalls, the pier surface and the deformed channel bottom were treated as no-slip boundaries. To account for the small unresolved roughness in the bed-surface mesh, the modified eddy viscosity at the wall was calculated such that the nondimensional bed roughness was close to 100. The boundary condition implementation for rough walls is described in Spalart (2000) and Zeng et al. (2008). The free surface was treated as a shear-free rigid lid, which is justified as the free surface deformations in the experiment were negligible $\left(\Delta z/H<0.04\right)$ . The maximum deformations $\left(\Delta z\sim 0.02\text{m}\right)$ in the experiment were recorded at the upstream stagnation point. The use of a rigid lid approximation accounts to a certain degree for the effect of the increase of the free surface level near the stagnation point on the flow, as the pressure distribution at the free surface is not uniform. The inferred free surface superelevation at the stagnation point was within 20% of the experimental value. The rigid lid approximation is widely used for simulating flows with ${\mathsf{F}}_H<0.5$ (e.g., see Paik et al. 2007). At the outflow, a convective boundary condition was used. The inflow conditions corresponded to fully developed turbulent channel flow and contained resolved turbulence fluctuations. The time step was $0.025D/U$ .

One issue is the capability of DES to accurately capture the dynamics of the unsteady coherent structures in the flow past surface-mounted bodies. The related DES-LES study of Koken and Constantinescu (2008a, 2009) for vertical wall abutments with incoming fully turbulent flow showed that the mean-flow and turbulence statistics predicted by a DES simulation at a channel Reynolds number of 18,000 were very close to LES results at the same Reynolds number, despite the different numerics in the two codes and, of course, use of a DES SGS model instead of the dynamic Smagorinsky model.

The largest values of $\overline{p^{\prime 2}}$ at the bed [Fig. 14(a)] are observed at the base of the pier $\left(r/D<0.7\right)$ . They are induced by (1) the strong temporal variations in the intensity of the downflow; (2) the movement of the legs of the necklace vortices; (3) the wrapping of the legs around the sides and the downstream face of the pier followed by detachment of patches of vorticity from the downstream part of the legs; and (4) the eddies shed inside the SSL in the vicinity of the bed. As the attached boundary layer on the pier separates at $\varphi \sim 180°$ , most of the SSL eddies in the near-bed region remain close to the symmetry plane. While some of these eddies are convected away from the pier, others are entrained in the region between the two SSLs and are pushed back toward the pier. The $\overline{p^{\prime 2}}$ levels inside most of the scour hole around the pier are also significantly higher (by at least two times) than those observed on the channel bed both upstream and downstream of the scoured region.

Fig. 14(b) shows the nondimensional bed-friction velocity, $u_{\tau }/U$ , distribution in the mean flow. As the governing equations were integrated through the viscous sublayer, $u_{\tau }$ was calculated using the definition ( $u_{\tau }={\left(\tau /\rho \right)}^{1/2}$ , where the bed shear stress is $\tau =\rho \left(\nu +{\nu }_t\right)U_1/\Delta n_1$ , $U_1$ is the magnitude of the velocity component in the direction parallel to the bed at the first point off the bed and $\Delta n_1$ is the wall-normal distance of the first point off the bed situated within the viscous sublayer). The bed-friction velocity is amplified beneath the region where the core of the MV and its legs oscillate. The largest values are observed close to the junction line between the pier and the bed for $50°<\left|\varphi \right|<160°$ The main reasons for the high values of $u_{\tau }$ close to the base of the pier are similar to those responsible for the amplification of $\overline{p^{\prime 2}}$ in the same region. Overall, the distribution of $u_{\tau }$ in the mean flow is similar to that predicted in case LR. By contrast, in the case of a flat bed, we observed that the maximum $u_{\tau }$ in the mean flow occurred at $\varphi \sim 45°$ on both sides of the pier. The amplification of $u_{\tau }$ in the region of strong flow acceleration is consistent with experiments that reported scour is initiated in this region. In the flat-bed case, $u_{\tau }$ is also amplified beneath the MV but the amplification is smaller than that observed in the acceleration region centered around $\varphi =\left|45°\right|$ .

Animations of the instantaneous distributions of the bed shear stress, $\tau$ , for the HR case show a high temporal variability of this quantity inside the scour hole. Though, at most times, the largest values of $\tau$ are observed in the regions where the bed shear stress is high in the mean flow, there are time instances when the largest values occur at a short distance from the back of the pier. Such an event is visualized in Fig. 15(a) (see arrow pointing toward the strong elongated patch of high $\tau$ values forming on the right side of the pier). Then, the streak of high bed shear stress starts being convected away from the pier while remaining on the side where it originated [Fig. 15(b)]. The main reason for the formation of these streaks of high bed shear stress is the movement of one of the legs of the necklace vortices toward the back of the pier followed by the detachment of a streak of vorticity from the leg (see discussion of Figs. 7 and 8). As this streak is convected away from the pier, its vorticity is high enough to produce a significant amplification of $\tau$ beneath it. Such a large amplification of $\tau$ behind the pier was not common in the LR case.

The dynamics of the large-scale coherent structures in the flow past a circular pier with a deep scour hole in a relatively shallow channel was investigated using DES for a flow condition typical of the upper range of Reynolds numbers at which flume studies of local scour around bridge piers are conducted under clearwater scour conditions. The flow mechanisms responsible for sediment entrainment and transport during the later stages of the scour process at a circular pier for ${\mathsf{R}}_D={10}^4-\left(2\times {10}^5\right)$ include bimodal, large-scale oscillations of the MV. Additionally, the curving of the legs of the necklace vortices toward the back of the pier followed by the detachment of a streak of vorticity and its convection away from the pier at a small distance from the bed (as depicted in Figs. 7 and 8) provide for a mechanism for sediment removal behind the pier. The topology of the flow in the near wake was similar in the low and high Reynolds number cases. In particular, the presence of strong upwelling motions near the symmetry plane behind the pier $\left(\varphi =180°\right)$ and the relative shallowness of the channel impeded the interaction between the SSLs and suppressed the formation and regular shedding of large-scale roller vortices in the near wake. Finally, in both cases the angle at which the attached boundary layers on the sides of the pier separate was found to increase from $90°$ in the upper part of the channel to close to $180°$ in the vicinity of the bed.

The DES and LES simulations revealed that some scale effects are present over the range of Reynolds numbers simulated. For example, differences were observed between the two cases in terms of: (1) the distributions of $k$ and $\overline{p^{\prime 2}}$ in the region where the large-scale bimodal oscillations of the MV were present; (2) the preferred mechanism associated with the transition between the two modes of oscillation; and (3) the probability the streak of vorticity detaching from the leg of the necklace vortex will be convected laterally or behind the pier. In particular, in the high Reynolds number case $k$ was strongly amplified in the region where patches of vorticity of opposite sign to that of the vorticity inside the core of the MV are ejected away from the bed, and the coherence of the streaks of vorticity detaching behind the pier was higher. In the high Reynolds number simulation the largest values of $\overline{p^{\prime 2}}$ in the HV region were observed in the region where the eddies convected with the downflow reached the scour hole. This means that at high Reynolds numbers the pressure fluctuations play a critical role in the entrainment of sediment near the junction line between the upstream face of the pier and the bed. Another important finding of this study is that the nondimensional values of $k$ and $\overline{p^{\prime 2}}$ in the regions of high turbulence amplification associated with the MV, the SSLs and the near wake, were consistently lower in the high-Reynolds-number simulation. This scale effect should be accounted for when results from experiments and numerical simulations at lower Reynolds numbers are used to calibrate scour prediction formulas. Scour prediction formulas calibrated using data mostly from low Reynolds number experiments tend to overpredict scour when used to predict scour at field conditions. The fact that the amplification of the turbulence in nondimensional terms is smaller at high Reynolds numbers suggests that Reynolds number effects are one of the reasons for the observed overprediction of scour in the field. Results of the present investigation do not allow a clear estimation of the importance of this effect compared to other factors responsible for the overprediction of scour in the field.

The present study considered a circular cylinder as pier shape. More complex shapes induce additional complexities of the flow field. The present simulation approach lends itself to investigating variations in flow structure and turbulence, as well as scour mechanisms with differences in such important variables as angle of attack, flow shallowness, pier width, and pier substructure. The presence of sharp edges (e.g., rectangular pier form), for instance, can change dramatically the angle at which the attached shear layers separate in the near-bed region (Kirkil and Constantinescu 2009). In turn, changes occur in the topology of the flow and the dynamics of the large-scale coherent structures in the near-wake region.

The writers thank the National Center for High Performance Computing in Taiwan and Dr. W. H. Tsai for providing the computational resources needed to perform some of the simulations as part of a collaboration program between the two institutions.