## 1 Introduction
The pier Reynolds number (\(R_{\rm{e}D}\)) is calculated based on the pier diameter (\(D\)) and the depth-averaged approach velocity (\(U_{o}\)). The \(DES\) was used for piers with a high pier Reynolds number, which lies outside the range of well-resolved \(LES\)[36; 37]. The \(DES\) is also used for rectangular piers with a high aspect ratio. The obtained results showed high bed shear stress values extended for a distance ranged between a 6 and 8 pier width upstream from the pier. Moreover, the nondimensional pressure root-mean-square (\(RMS\)) fluctuations at the bed were one order of magnitude higher than the simulated ones for circular piers [38]. The LES captured bimodal random oscillations in the necklace vortices forming the horseshoe system during all scour stages [36; 37]. It was concluded that the flow remains subcritical near the bed surface due to the reduction in the boundary layer, for both a high and low pier Reynolds number [39]. The sub-particle-scale (SPS) turbulence closure scheme was used in the three-dimensional Incompressible Smoothed Particle Hydrodynamics (ISPH) erosion model to account for turbulence scales that are below the sediment particle size. The three-dimensional ISPH model was utilized to simulate the local scour at a relatively large cylindrical pier by using the turbidity water particle concept (TWP). The obtained scour hole depths had a suitable degree of accuracy for engineering purposes in some cases, and significant errors in other cases. It was concluded that obtaining a highly reliable scouring effect is very difficult even with experimental studies [40]. The eddy viscosity and turbulent shear stress were estimated by applying the SPS turbulence model with the Smagorinsky approach. The accurate investigation of the SPS turbulence closure scheme especially at the interfacial boundaries was not possible due to the lack of detailed turbulence data [41]. The continuous improvement in turbulence modeling, especially within the roughness layer in the macroscopic simulation of rough boundaries, will enhance the _CFD_ modeling for the bed shear stress and accordingly, the sediment initiation of motion. A proposed three layer mixing length model was developed to represent the eddy structure sizes in the porous bed, the roughness layer, and the free flow zones in the case of channel flows over and within natural porous beds [42]. In short, _CFD_ techniques are beneficial tools to understand the complexity of flow fields around bridge piers. Despite some successful scour modeling attempts, _CFD_ techniques are characterized by their limited ability to simulate scour and erosion processes; as a result, their fidelity in scour depth estimation is limited [43]. The two main objectives of the current research are:  
* Selecting an appropriate turbulence model with personal computer computational resources that can simulate the case in hand with acceptable accuracy, since no agreement between researchers has been achieved yet;
* Tracking the flow structure and bed shear stress variations during the development of a scour hole with the maximum achievable number of development stages. The morphological dimensions for the development stages (seven stages) of the scour hole were collected from previous studies.