Stationary single waves in turbulent free‐surface flows

Stationary near‐critical turbulent open‐channel flow is considered. Allowing for locally enlarged roughness of the inclined bottom leads to the development of a stationary single wave. The waveform of the free surface is the key interest of this work. The flow is assumed to be fully developed far upstream and far downstream. In order to obtain the shape of the surface elevation as a result of the full equations of motion, an iteration procedure is used, which was successfully applied to the related problem of the undular hydraulic jump [1]. In addition to the classical solitary wave, the so‐called single wave of the second kind is studied numerically. The numerical results are compared with theoretical predictions and experimental data.


Introduction and problem description
Ensemble-averaged quantities are denoted by overbars. The volumetric mean flow velocity and the surface elevation in the fully developed reference state are denoted byūr andhr, respectively. The acceleration due to gravity is indicated by g. The x-axis is in the bottom plane, while the y-axis points upwards. The gray box indicates the region of enlarged bottom roughness of length , leading to a piecewise constant distribution of the friction coefficient, with c f > c f r . This phenomenon has been studied repeatedly in the past, performing asymptotic analyses as well as experiments [2][3][4]. According to analytical investigations [4] two types of waves can occur under the before mentioned circumstances, namely a wave that resembles the well known solitary wave, which will be denoted as single wave of first kind, and the so-called single wave of second kind. The shapes of both kinds of single waves differ from each other not only in terms of position and amplitude of the wave crest but also in terms of the decay behavior upstream as well as downstream.

Iteration procedure and results
Due to the near-critical flow conditions it may be difficult to satisfy the boundary conditions at the free surface with finitevolume or finite-element methods. In order to circumvent these difficulties, it is proposed to apply an iteration procedure that is based on an asymptotic expansion for upstream Froude numbers Fr :=ū r / gh r approaching the critical value 1, i.e.
referred to the hydrostatic pressure in the reference state with constant density ρ, is prescribed at the surface, serving as a basis for the iteration procedure. The Reynolds number, defined in terms of the friction velocityū τ r = gαh r for the fully developed flow, is assumed to be very large, i.e. Re τ := gαh 3 r /ν → ∞, where ν is the constant kinematic viscosity. An asymptotic analysis of the governing equations in non-dimensional form, i.e. continuity equation and Reynolds-averaged Navier-Stokes equations, in terms of ε, yields the following steady state version of an extended Korteweg-de Vries (KdV) equation for the dimensionless surface,h(x)/h r = 1 + εH 1 (X) + O(ε 2 ): In eq. (1) primes indicate derivatives with respect to the contracted non-dimensional coordinate X = 3 √ εx/h r . The parameter β := (1/3)αε −3/2 characterizes the effect of turbulent dissipation, and Γ(X) = (c f − c f r )/(3εc f r ) = O(1) is a forcing function due to the enlarged bottom roughness. As described in [1], the iteration procedure starts with numerically solving eq. (1) for P S ≡ 0 in order to obtain the asymptotic 2 of 2 1 as a first guess of the free surface elevation. In the next step, the full equations of motion are solved numerically using the commercial software FLUENT with the standard k-ε turbulence model. Using a more elaborate turbulence model has no significant effect on the solution, cf. [1]. The computed pressure distribution P S (X) at the prescribed surface elevation H (1) 1 differs from the value 0 that is required by the dynamic boundary condition. To determine a corrected free-surface shape H (2) 1 , eq. (1) is solved with P S (X) from the preceding FLUENT computation. The function H (2) 1 is then used to determine the solution of the full equations of motion in the next iteration step. The procedure is repeated until the bounds of accuracy for the surface pressure and the surface elevation, respectively, are reached. For the present case the following bounds were considered as sufficient: max{|P S (X)|} < 10 −3 and max{|H (n+1) 1 − H (n) 1 |} < 10 −4 . This method was applied to the solutions of the first and the second kind, respectively, for parameter values matching the experimental values according to [2].
In Fig. 2 the numerical solution of the full equations of motion (red curve) is compared with the asymptotic solution (blue curve) and with experimental data [2]. After 28 iteration steps the convergence criteria were satisfied. The position of the wave crest is well predicted by both the asymptotic and the numerical solution. The numerically obtained amplitude shows a higher discrepancy from the measured surface amplitude, but the numerical result is in better agreement with the experimental data in terms of wave width and shape. Presently, the reason for that discrepancy is not clear. As can be seen in Fig. 2, a negative wave tail was measured at a certain position behind the region of enlarged bottom roughness, leading to a weak downstream decay, like in the numerical solution.  Fig. 3: Non-dimensional first-order surface elevation of the second kind for the same set of parameters as in Fig. 2; blue: asymptotic solution [4], red and green: non-convergent numerical solutions after 15 and 28 iteration steps, respectively.
In Fig. 3 the starting point of the iteration routine was the solution of the second kind of eq. (1) with P S ≡ 0, characterized by a long shallow tail and a relatively small amplitude (blue curve). A conspicuous growth of the wave height can be observed with an increasing number of iteration steps, cf. red and green curves. As a consequence of a non-decreasing pressure perturbation at the surface during the iteration procedure, the wave continues to grow, until the assumption P S = O(ε 2 ) 1 is no longer valid. As a result, eq. (1) is no longer solvable after the 28 th iteration step. At the same time, the wave crests are moving upstream. The surface elevation seems to approach a solution of the first kind, which is, however, unstable [2]. Certainly the problem deserves further investigations, in particular concerning the stability of the solution of the second kind. In that context it may be of interest to note that the present problem is mathematically equivalent to the problem of near-critical flow over a bump of small height [5]. In the latter case, however, surface elevations in accord with solutions of the second kind have indeed been found both theoretically and in experiments [5].