Near-critical turbulent open-channel ﬂows over bumps and ramps

Steady two-dimensional turbulent free-surface ﬂow in a channel with a slightly uneven bottom is considered. The shape of the unevenness of the bottom can be in the form of a bump or a ramp of very small height. The slope of the channel bottom is assumed to be small, and the bottom roughness is assumed to be constant. Asymptotic expansions for very large Reynolds numbers and Froude numbers close to the critical value Fr = 1, respectively, are performed. The relative order of magnitude of two small parameters, i.e. the bottom slope and ( Fr − 1 ) , is deﬁned such that no turbulence modelling is required. The result is a steady-state version of an extended Korteweg–de Vries equation for the surface elevation. Other ﬂow quantities, such as pressure, ﬂow velocity components, and bottom shear stress, are expressed in terms of the surface elevation. An exact solution describing stationary solitary waves of the classical shape is obtained for a bottom of a particular shape. For more general shapes of ramps and bumps, stationary solitary waves of the classical shape are also obtained as a ﬁrst approximation in the limit of small, but nonzero, dissipation. With the exception of an eigensolution for a ramp, an outer region has to be introduced. The outer solution describes a ’tail’ that is attached to the stationary solitary wave. In addition to the solutions of the solitary-wave type, solutions of smaller amplitudes are obtained both numerically and analytically. Experiments in a water channel conﬁrm the existence of both types of stationary single waves.

become singular as the critical Froude number 1 is approached. We shall further discuss the point in Sect. 4 in the context of an extended KdV equation. Furthermore, no upstream influence of the obstacle is predicted by the classical approach, which is in contrast to observations on near-critical flows, including the experiments to be presented in this paper.
With the exception of the stability analysis given in Appendix D, the present analysis concerns steady flow. The same is true for the following discussion of the pertinent literature, leaving propagating waves mainly unnoticed.
For inviscid flow, a conventional approach is to solve the equations of potential flow or the Euler equations, subject to boundary conditions at the bottom and at the free surface. A comprehensive survey on both classical work and more recent research is given in [11], Ch. 3. As examples that are of some relevance to the present work, we mention the following ones. Forbes and Schwartz [20] considered a semi-circular obstacle, analysed linearized equations, and, in addition, provided numerical solutions of the nonlinear problem. Numerical solutions for a semi-circular obstacle were also investigated by Vanden-Broeck [67], and a non-uniqueness was found, with one solution corresponding to a perturbation of the uniform stream, the other one to a perturbation of a solitary wave. Numerical solutions for a triangular obstacle were given by Dias and Vanden-Broeck [16]. Rapid changes in the channel bottom give rise to particular problems, of course. An example is a step in the channel bottom; see [2].
Under certain assumptions, the full equations of motion of inviscid flow can be reduced to the Korteweg-de Vries (KdV) equation for the free surface, or to various modifications and extensions of that famous equation, cf. the monographs [51,63] as well as the survey articles [15,26,29]. In their work mentioned already above, Binder et al. [2] solved the full nonlinear problem numerically, and, in addition, provided solutions of a forced KdV equation for the weakly nonlinear problem. The semi-circular obstacle found again particular interest, see [17], a non-uniqueness was observed, with two symmetric solitary waves as solutions, and a solution with a cnoidal wave downstream for upstream supercritical flow was found. Nonlinear wave trains (resembling unsteady undular bores) upstream and downstream of a locally steady surface elevation over the obstacle were obtained in [27,28]. The interesting phenomenon of "shelves" was described in [39], whereas steady or oscillatory "tails", attached to solitary waves, were predicted in [37] for higher-order KdV equations. Nonuniqueness was observed in [1,22], while the stability of the solutions was investigated in [4,5,13], among others. Finally, it might be of interest to mention that coupled KdV equations describing resonant flow over an uneven bottom have also been studied [44]. Concerning the effects of viscosity it may be mentioned that weakly dissipative free-surface flows using the classical potential-flow approach, based on linearized Navier-Stokes equations, were investigated in [18,19]. Remarkably, non-local viscous terms in the boundary condition at the bottom were found in [19].
Turbulent flow, which is the rule rather than the exception in open-channel hydraulics, is, of course, associated with dissipation. A variety of dissipation models can be found in the literature; cf. [11], Ch. 5, for a survey that gives also many important details. Several authors applied the so-called Chezy model to account for bottom friction, i.e. the wall shear stress is assumed to be proportional to the mean-velocity square; cf. [45,46]. In the present context, it ought to be mentioned that Grillhofer [24] added a friction term of the Chezy type to the basic equations for plane steady inviscid flow, expanded the equations for Froude numbers close to 1 and obtained an extended KdV equation for the surface elevation with a linear dissipation term. Note, however, that linear terms in extended KdV equations to describe bottom friction and dissipation are sometimes attributed to the "Rayleigh model of dissipation", cf. [55]. In view of results to be presented in the present work, it may be of interest that the possible appearance of tails in case of small, but non-vanishing dissipation was discussed in [56]. It ought to be mentioned, however, that other dissipation models have been proposed as alternatives, e.g. in [6][7][8][9][10]. The model offered in [8] also accounts for the effects of streamline curvature.
There is, of course, always the possibility to solve the full equations of motion for turbulent flow. Presently available commercial codes, e.g. FLUENT or FLOW-3D, allow the computation of turbulent shear flows with free surfaces based on the Reynolds-averaged Navier-Stokes equations with various turbulence models. However, finding the solution of a particular problem is not always an easy task when the Froude number is close to the critical value 1 and the flow to be described contains both supercritical and subcritical regions; see [32]. To cope with those difficulties, an iteration method based on an asymptotic analysis of near-critical turbulent free-surface flow was proposed in [61] and applied to the undular-jump problem. Another example of a successful solution of the Reynolds-averaged Navier-Stokes equations for undular jumps was recently given in [58], using the volume of fluid method. It is remarkable that a one-dimensional flow model based on the depth-averaged Navier-Stokes equations with Reynolds stress terms can give results that are comparable with the solutions of the two-dimensional flow problem, see [11], Section 5.3.5. At a higher level of computational effort, large-eddy simulation (e.g. [31]) and even direct numerical simulations (e.g. [30,40,50]) of turbulent free-surface flows are also available. From the view of the present work, the DNS data due to [30,40] will be of importance as they provide justification for the assumption of vanishing Reynolds shear stress at the free surface, cf. below.
The aim of the present work is to obtain asymptotic solutions for near-critical turbulent open-channel flow over bumps and ramps without applying turbulence models. An extended KdV equation, containing a linear dissipation term and a forcing term, was obtained in [60] for a channel bottom with a region of enlarged roughness, providing the possibility to predict, apparently for the first time, stationary solitary waves in turbulent open-channel flow. In the present work, the analysis due to [60] is modified to account for an uneven bottom, e.g. for bumps or ramps, while keeping the bottom roughness constant. The main result of the analysis is an equation for the averaged surface height as a function of the longitudinal coordinate, but other flow quantities of interest are expressed in terms of the surface elevation.
For the previous case of a plane bottom with a region of enlarged roughness, the stationary solitary waves were observed in laboratory experiments and reasonable agreement between predictions and measurements was found [62]. Those favourable comparisons made it appear feasible to verify the theoretical predictions by experiments also for bumps. However, the experiments led to surprising observations that will be described in Sect. 7 below.

Reference quantities and non-dimensional parameters
Steady two-dimensional turbulent open-channel flow is considered (Fig. 1). Apart from small perturbations of characteristic length l and elevation b, which varies in longitudinal direction, the channel bottom is plane with small constant slope α. The order of magnitude of b will be fixed in course of the analysis. For now, it suffices to assume that b is much smaller than the depth of the liquid in the channel. A Cartesian coordinate system is chosen such that the x-axis is in the bottom plane, while the y-axis points upwards. The flow velocity components in the x,y coordinate system are u and v, respectively. Ensemble-averaged quantities will be denoted by an overbar; e.g.h stands for the ensemble-averaged surface height. Since steady flow is considered, the ensemble-averaged quantities are in accord with time-averaged ones. Primes will indicate fluctuations around the average.
It is assumed that the flow would be fully developed 1 if the bottom were exactly plane. For bumps or ramps, this implies that the flow is indeed fully developed far upstream and far downstream. Note that in case of a ramp the surface elevation is slightly larger far downstream than far upstream, cf. Fig. 1. However, the difference will turn out to be small of higher order in the asymptotic analysis.
For introducing non-dimensional variables, the fully developed flow over a plane bottom with slope α is chosen as a reference state, which is denoted by the subscript r. Thus, the averaged surface height of the fully developed flow far upstream,h r , serves as reference length, and the y-coordinate will be referred to it. The x-coordinate, however, will be contracted with a small parameter δ, which is chosen in a way that is appropriate for describing solitary waves of small amplitudes. From [60], it is known that a suitable choice is where the non-dimensional parameter ε characterizes the amplitude of the wave, while the coefficient 3 has been introduced in the interest of simplifying the coefficients of the final equations. A definition of ε in terms of the Froude number will be given below, cf. (17). The non-dimensional Cartesian coordinates X, Y and the non-dimensional bottom elevation B are then defined as follows: For the purpose of serving as reference velocity, the volumetric mean velocityū r is introduced, i.e.ū r = V /h r , whereV denotes the volume flow rate per unit width of the channel. This choice of the reference velocity is of importance for keeping the analysis free of turbulence modelling, see [25]. Based on that reference velocity, the Froude number is defined as where g is the acceleration of gravity. Concerning further dimensionless quantities, the pressurep is referred to the averaged hydrostatic pressure at the bottom of the channel in the fully developed flow, i.e. gρh r , with ρ denoting the constant density of the fluid. The Reynolds stresses are referred to the averaged bottom shear stress in the fully developed flow, which may be written as ρu 2 τr , where u τr is the reference friction velocity. Since the wall shear stress balances the tangential component of the gravity force in fully developed flow, u τr is given by As the bottom slope, α is assumed to be very small, cf. above, sin α was replaced by α in Eq. (4). Using the reference quantities defined above, non-dimensional variables are introduced as follows: The effects of surface tension will be neglected. Concerning viscosity effects, very large Reynolds numbers will be considered. Following common practice in the asymptotic analysis of turbulent flow, the Reynolds number is defined in terms of the reference friction velocity, i.e.
where ν is the kinematic viscosity of the liquid.

Governing equations for large Reynolds numbers
For very large Reynolds numbers, the flow field is composed of two layers, i.e. the defect layer and the thin viscous wall layer at the bottom. Concerning the latter, a universal solution is known to exist for steady flow, cf. [21], or [59]. Thus, it is sufficient to consider only the defect layer in what follows.

Equations of motion
In terms of the present non-dimensional variables, the continuity equation for the mean flow reads Here, and in what follows, derivatives are indicated by subscripts, whenever it is convenient.
Since viscous stresses are negligible in the defect layer, the momentum equations for the mean flow become Here, and in what follows, sin α and cos α are replaced by α and 1, respectively.

Boundary and matching conditions
With regard to the boundary condition at the bottom, it is assumed that the thickness of the viscous wall layer is much smaller than the maximum bottom elevation, b max , i.e. ν/u τr b max or For a given bottom elevation, Eq. (9) is a condition for the Reynolds number, which will be specified more precisely as soon as the order of magnitude of B max has been fixed, see Sect. 4 below. If Eq. (9) is satisfied, the thin viscous wall layer will be attached to the bottom elevation. Otherwise, the bottom elevation would be submerged partially, or even totally, in the viscous wall layer, requiring a special treatment; cf. [38,66].
With the viscous wall layer attached to the bottom elevation, the appropriate boundary condition for the defect layer is the classical tangential-flow condition, i.e.
In order to match the defect-layer solution to the wall-layer solution, the logarithmic overlap law is applied, following [59], p. 536. The following expression for the non-dimensional surface velocity is obtained: where U τ andC(X ) are defined by the relations andC respectively. Furthermore, κ is v. Kármán's constant, while the value of C + depends on the non-dimensional roughness, k + s , with where k s is the sand roughness height, see [59, pp. 528-530]. In the present analysis, the sand roughness is assumed to be constant over the bottom. Thus, C + is also constant. Neither the value of C + nor the value of v. Kármán's constant matter for what follows.
It should be noted that (11) does not imply the assumption of a logarithmic, or any other particular, velocity profile across the defect layer. (11) is just a result of applying the matching condition to the defect layer, which is characterized by very small perturbations of the volumetric mean velocity.
Although it is not of relevance to the further analysis, it may be of interest to note that U τ is exactly equal to the non-dimensional friction velocity only in case of a plane bottom, i.e. B(X ) ≡ 0; otherwise, U τ differs slightly from the local non-dimensional friction velocity, i.e. by a higher-order term in the asymptotic expansion presented in Sect. 4.
Since there is a wall-pressure gradient at the bottom due to the X-dependence of the water depth, it is in order to check whether the logarithmic overlap law is applicable. An appropriate analysis is given in "Appendix A".
At the free surface, kinematic and dynamic boundary conditions have to be satisfied. The kinematic boundary condition is prescribed in the conventional form, i.e.
which may be interpreted as a definition of the ensemble-averaged surface. With respect to the dynamic boundary conditions, we also follow the conventional approach, see [57]. Since surface tension is neglected and viscous stresses are negligibly small in the defect layer, only the Reynolds stresses have to be taken into account. Thus, the required continuity of the total stresses at the free surface leads to the following dynamic boundary conditions with respect to the X-and Y -directions, respectively: where ϑ is the inclination angle of the ensemble-averaged free surface with respect to the horizontal, i.e. tan ϑ = δH X . The application of "conventional" boundary conditions at the free surface rests on the assumption of vanishing Reynolds shear stress at the free surface. This assumption has been subject to some discussion; see [60] and the references given there. It may well be that the apparently unsettled question of zero versus nonzero Reynolds shear stress at the free surface is associated with problems concerning the ensemble-averaged kinematic boundary condition at free surfaces, as noted in [65] cf. also [43]. Furthermore, recent investigations [42] show a three-dimensional character of the turbulent flow at free surfaces; however, implications with regard to boundary conditions for the bulk flow seem to be uncertain at present. For the present analysis, it may be seen as reassuring that the application of the conventional boundary conditions (16a) and (16b) in previous work has led to reasonable agreement with measurements of surface elevation as well as shear stress distribution for both undular jumps [25,35,36,61] and stationary solitary waves [62].

Asymptotic expansion for near-critical flow
Since solitary waves with small amplitudes move in an inviscid quiescent fluid with a velocity that corresponds to a slightly supercritical Froude number, it is natural to expand the governing equations in terms of a small parameter ε that is defined by the relation Here, and in what follows, coefficients like 3/2 are introduced in order to obtain the results in a convenient form.
In general, the small parameter ε is independent of the slope α that has been introduced above as a small parameter. However, it was shown in earlier work [25,60] that the analysis can be kept free of turbulence modelling if α and ε 2 are assumed to be of the same order of magnitude. The basic idea is to fix the relative size of the two small parameters ε and α such that, on the one hand, the leading terms due to turbulence, though small, are retained in the analysis, while, on the other hand, the magnitude of the Reynolds stresses is sufficiently small to keep the analysis free of turbulence modelling. Thus, following [25,60], the coupling parameter A is introduced with the relation Note that the fully developed flow with small bottom slope α requires a particular, very large value of the Reynolds number. Thus, (18) can also be seen as a coupling condition for the nearly critical Froude number and the large Reynolds number. A further important assumption concerns the height of the bottom elevation (e.g. a bump or a ramp). It is chosen such that the terms due to the bottom elevation will become of the same order of magnitude as the terms due to dissipation. This is accomplished by introducing the parameter β with the definition and assuming that the non-dimensional bump or ramp height defined in (2) is of the order of magnitude B = O(βε 2 ) = O(ε 5/2 ). Thus, the non-dimensional height of the bottom elevation is written in the following form: For the further analysis it is of interest that the order of magnitude of B according to (20) allows an analysis that remains free of turbulence modelling. It should also be noted that substituting (20) and (19) into Eq. (9) gives α √ ε Re τ 1 as the condition for sufficiently thin viscous wall layers at the uneven bottom.
That condition, or an equivalent one, has to be observed when experimental verification of the theoretical predictions is sought; cf. Sect. 7. For bumps, the bottom elevation vanishes far upstream and far downstream. In case of ramps, the bottom elevation vanishes only far upstream, whereas it attains a nonzero constant value far downstream. This leads to the following conditions for ψ: with ψ ∞ = 0 for a bump, and ψ ∞ > 0 for a ramp. The dependent variables are now expanded in terms of powers of ε as follows: and analogous expansions for the velocity components, the pressure and the Reynolds stresses. The expansions are then introduced into the governing equations given in the preceding section. The procedure is a bit cumbersome and not always straightforward, as the averaged Navier-Stokes equations remain incomplete without a closure model. The details of the analysis, which follows the analogous analysis for a plane bottom with a region of enlarged bottom roughness [60], can be found in "Appendix B". For here it may suffice to mention that the main difference between the present asymptotic expansion and the previous one stems from the boundary condition (10) together with (20), which, of course, has no equivalent in case of a plane bottom. Eventually, the following first-order results are obtained: where H 1 has to satisfy the ordinary differential equation The term U (Y ) appearing in (23) is the non-dimensional "velocity defect", which is defined such that the local value ofŪ in the fully developed flow far upstream differs from its volumetric mean, i.e. 1, by the amount √ α U . This quantity can be considered as known from experiments, e.g. [52], but, in the present analysis, it does not affect any other variable exceptŪ . Note that the distribution of the first-order perturbation of the Reynolds shear stress over the channel cross section cannot be obtained in the framework of the present analysis; only the values at the bottom (Y = 0) and at the undisturbed surface level (Y = 1) can be determined with (26) and (27), respectively. According to Eqs. (23) to (27), all first-order perturbation quantities are expressed in terms of the first-order surface elevation, H 1 (X ). The latter is to be determined by solving (28), which can be recognized as a steadystate version of an extended KdV equation. The extension consists of two parts, i.e. the linear term β H 1 , which characterizes dissipation, and the "forcing" term −βψ X , which is a given function of X, making the equation an inhomogeneous one. If the term β H 1 is dropped, (28) describes inviscid flow. As already anticipated above, U does not appear in (28), so that knowledge of the velocity defect is not required. For later use it will be of interest to observe that replacing ψ X in (28) by Γ leads to the extended KdV equation derived previously for a plane bottom with a roughness that varies along the bottom [60,62]. However, the formal similarity between the results should not obscure the fact that the forcing terms differ in origin. Γ results from an expansion of the logarithmic overlap law (11) with C + (X ) according to the roughness as a function of X , whereas ψ X follows from expanding the boundary condition (10) with (20) for the shape and height of the bump or ramp.
It may be worth mentioning that dropping the term H 1X X X in (28) leads to an equation that is equivalent to the one-dimensional flow approximation of classical hydraulics [24]. 2 Obviously, the remaining equation is singular as the surface elevation tends to the critical value H 1 = 1. This shows the inability of the onedimensional flow approximation of classical hydraulics to describe the near-critical flow considered in the present paper.
Since we consider flows that are fully developed both far upstream and far downstream, the solutions of (28) have to satisfy the boundary conditions The boundary condition of vanishing surface elevation as X → + ∞ is trivial for a bump, as the bottom elevation vanishes far downstream in this case. For a ramp, however, there remains a nonzero free-surface elevation that is of the order of the bottom elevation as the fully developed flow far downstream is attained.
With the bottom elevation being as small as defined by (20), the free-surface elevation far downstream is small of higher order and does not affect the first-order boundary condition (29). Integrating (28) from −∞ to +∞ and accounting for the boundary conditions (29), the following integral relation is obtained: Equation (30) represents the conservation of momentum for steady flow, see "Appendix C". By integrating (28) after multiplication with H 1 and accounting again for the boundary conditions (29), one obtains another integral relation that reads as follows: If the bottom is uneven only in a region of finite length (in terms of the contracted coordinate X ), the solutions of (28) subject to the boundary condition (29) decay exponentially, i.e. H 1 ∼ exp(k X), as X → ± ∞. Since (28) is of third order, k has to satisfy a cubic equation that reads As at least two real roots of (32) are required for describing a solution that decays both upstream and downstream, a necessary condition for the existence of such a solution of (28) is β ≤ 2/3 √ 3 ≈ 0.385, i.e. dissipation must not be too strong. For β < 2/3 √ 3 the cubic equation (32) has two negative real roots, k 1 < k 2 < 0, and one positive real root, k 3 > 0. For very small values of β, (32) has the roots k 1 = −1+O(β); k 2 = −β +O(β 3 ); Thus there are two types of solutions that differ from each other with respect to the decay downstream, i.e. as X → ∞. One type of solution is characterized by a strong decay according to the root k 1 , the other one by a weaker decay according to k 2 . Examples will be given below.

The stationary solitary wave as an exact solution
For a ramp of the shape an exact solution of (28) with boundary conditions (29) is i.e. the classical solitary-wave solution of inviscid flow, see Fig. 2. The term in parenthesis on the right-hand side of (28) is identically zero in this particular case. Thus, the solution (34) has the remarkable property of being valid for any value of the dissipation parameter β. In case of inviscid flow, i.e. β = 0, the right-hand side of (28) would vanish and the position of the stationary solitary wave would remain undetermined. In other words, it is the dissipation that fixes the position of the stationary solitary wave relative to the position of the ramp.  (34) in turbulent open-channel flow over a ramp of shape (33) Note that, in the present formulation, H 1 represents the elevation of the free surface with respect to the channel bottom including the bottom elevation due to the bump or the ramp, whereas the solitary wave in inviscid flow would require a plane bottom. However, as the height of the bump or ramp is assumed to be much smaller than the surface elevation, cf. (20), the distinction is not of relevance in the present case.

Solutions for small dissipation parameters
The dissipation parameter β is quite small in many practical applications as well as in pertinent laboratory experiments. For instance, comparisons with measurements were made for stationary undular jumps with β = 0.086 . . . 0.225 [64], β = 0.0563 [36], β = 0.00764 and β = 0.086 [35], among others, and for a stationary solitary wave with β = 0.127 [62]. Thus it appears appropriate to investigate the behaviour of the solutions for β 1.

Solutions of the first kind (solitary-wave type): general relations
In the limit β → 0, the first-order non-trivial solution of (28) with boundary conditions (29) is where it fixes the position of the maximum surface elevation in first order. The parameter X (0) m depends on the shape of the bump or ramp and can be determined from the integral relation (31). However, before that is done, it is advisable to check whether the integral relation (30) is satisfied. If the first-order solution (35) is inserted into the integral relation (30), one obtains ψ ∞ = 12. Thus, (35) is a uniformly valid solution of the extended KdV equation (28) with the homogeneous boundary conditions (29) only in case of ramps with the non-dimensional height ψ ∞ = 12. This value plays the role of an eigenvalue to the problem.
For bumps, as well as for ramps with ψ ∞ = 12, the integral relation (30) requires that a long, but shallow, tail is added to the solitary-wave solution (35). Then, (35) serves as a leading-order inner solution. An asymptotic analysis for β 1, analogous to the case of a plane bottom with a region of enlarged roughness considered in [62, Appendix II], gives the following leading-order solution for the tail (outer solution) in terms of the outer coordinate β X > 0: According to (36), the shape of the bump or ramp does not affect the tail. The tail depends, in leading order, only on the non-dimensional bottom elevation ψ ∞ far downstream, which is zero for a bump and nonzero for a ramp. If ψ ∞ > 12, the surface elevation is positive in the whole stationary wave, including the tail. If ψ ∞ < 12, the tail is a shallow trough. For ψ ∞ = 12 the tail vanishes, which is in accord with the eigensolution given above. If the dissipation term β H 1 is dropped in the extended KdV equation (28), the integral condition (30) is trivially satisfied, in first order, by the classical solitary-wave solution (35) in case of a bump of any shape, whereas it cannot be satisfied in case of a ramp. In the former case, the integral condition (31) gives X (0) m = 0 for a symmetric bump, i.e. the wave crest is exactly above the top of the bump in this particular case of inviscid flow. However, it should be noted that the solution for the inviscid case may not exist, though the solution of the full equation, i.e. including the dissipation term, does. A well-known related problem is the undular jump, which cannot be described as an inviscid flow, no matter how weak dissipation may be [25].

Solution of the first kind (solitary-wave type) for a plane ramp
An example that may be of practical interest is a plane ramp, i.e. a ramp with constant gradient ψ X = ψ ∞ /L in the region 0 < X < L. Note that the geometric height and length, respectively, of the ramp are b ∞ = ψ ∞hr α √ ε and l = Lh r /3 √ ε. The gradient of the ramp is b ∞ /l = 3ψ ∞ αε/L, which is as small as ε 3 . The integral relation (31) then gives the following expression for the position of the crest of the solitary wave: which is analogous to equation (5) of [48] for the problem of a plane bottom with enlarged constant roughness in the region 0 < X < L. A graph of (37) is shown in Fig. 3a. Real values of X Equation (38) defines a minimum non-dimensional height ψ ∞ for a given non-dimensional length L of the ramp, or, vice versa, a maximum length for a given height. If ψ ∞ = 12, which leads, in first order, to the stationary solitary wave without a tail, Eqs. (37) and (38) reduce to and respectively. Formally the same result was obtained for a plane channel bottom with enlarged constant roughness in a region of length L [62]. It was shown there in Appendix I that only the solution with the lower sign in (39) corresponds to a stable stationary wave. By analogy, it can be assumed that the same is true for (37). According to (40) stationary solitary waves without a tail exist for plane ramps only if the ramp length does not exceed a maximum value L max , which is given by L max ≈ 5.151 in first order. For constructing a uniformly valid solution including the non-vanishing tail according to (36) the secondorder inner solution is required. In general, the second-order analysis will be rather cumbersome. However, in the special case of a plane ramp, the results can be taken from [47], where the mathematically equivalent problem of a plane bottom with a region of constant enlarged roughness is considered. The parameter Γ L with the eigenvalue λ of [47] corresponds to ψ ∞ /L in the present problem. This gives the following second-order eigenvalue for the ramp height: m , according to (37) and (42), respectively. a Plane ramp with height ψ ∞ and length L. First-order eigenvalue ψ ∞ = 12. b Bump with isosceles triangular cross section of height ψ(0) and half-length L The second-order eigenvalue may be of importance for predicting the parameters for experiments with good accuracy. Since neither the second-order inner solution nor the uniformly valid solution will be used in what follows, we refrain from reproducing the rather lengthy results here, and refer to [47].

Solution of the first kind (solitary-wave type) for a bump with triangular cross section
Another case of interest is a bump with a cross section of the shape of an isosceles triangle in the region −L < X < +L. The piecewise constant gradient of the bump, ψ X , is then related to the non-dimensional height of the bump, ψ(0), by the equation ψ X = ±ψ(0)/L, with the upper and lower signs for −L < X < 0 and 0 < X < +L, respectively. The geometric height of the top of the bump is b t = ψ(0)h r α √ ε, while the geometric length of the bump is l = 2Lh r /3 √ ε. A similar bump was investigated in [16] based on potentialflow theory. Triangular bumps are also relatively easy to produce for the purpose of laboratory experiments; cf. below.
For turbulent flow, the integral relation (31) gives the following implicit relation for the position of the wave crest, It follows from (42) that a minimum height of the bump is required for the stationary solitary wave to exist, i.e.
Provided ψ(0) > ψ(0) min , (42) has two real solutions for X  [63], predicts stability. However, the experiments that will be described below indicate that the realization of the "stable" solution of the solitary-wave type is not an easy task and may even be impossible.

Solutions of the second kind
The solitary-wave solution (35) is associated with a surface elevation of the order of ε. However, (35) is not the only non-trivial solution of the extended KdV equation (28) with boundary conditions (29). The previously investigated case of a plane bottom with a region of enlarged roughness may serve as a guidance. Instigated by numerical solutions of the unsteady, transient flow problem [34], it was shown [47] that there are also solutions of (28), (29) with the leading term being of the order of β, i.e. describing free-surface elevations of the order of βε = O(ε 3/2 ). Those solutions describe waves that will be addressed as "stationary single waves of the second kind" in what follows. They are characterized, among others, by remaining supercritical in the whole flow field if β 1. Since the problem of a plane bottom with constant enlarged bottom roughness in a region of non-dimensional length L is mathematically equivalent to the problem of a plane ramp as defined in Sect. 5.2.2, first paragraph, the uniformly valid solution given in [47] may be re-written in terms of the present variables as follows: with The superscript in indicates that it concerns the contribution of the inner solution.
It is remarkable that the solutions of the second kind for a ramp, in contrast to the solutions of the first kind (solitary-wave type), exist only for non-vanishing dissipation, display weak downstream decay, and the height of the tail is of the same order of magnitude as the main part of the stationary wave.
For a bump, the behaviour of the solution of the second kind differs considerably from that given above for a ramp. According to [49], the inner solution displays strong decay and, therefore, is a uniformly valid solution of first order. To find a uniformly valid solution that comprises the shallow tail would require to determine the second-order inner solution and match it with the outer solution. However, the extended KdV equation (28) was derived by an expansion for small values of ε, and it seems doubtful whether it is still valid for terms being as small as the second-order terms in the inner expansion. Thus, we refrain from providing a uniformly valid solution of the second order and rather use numerical solutions for the comparison with experimental data.

Method of numerical solution
Even for simple geometries of the ramp or the bump, e.g. a plane ramp or a bump with an isosceles triangular cross section, the problem is governed by three independent parameters, i.e. the non-dimensional length L, the non-dimensional height ψ ∞ or ψ(0), and the dissipation parameter β. In the lack of a general theory, it would require an enormous computational effort to determine the regime of existence of solutions in the 3D parameter space. That was beyond the scope of the present work. However, numerical solutions have been obtained, if possible, for the purpose of comparisons with analytical solutions and experimental data, respectively.
Though the steady-state version (28) of the extended KdV equation is an ordinary differential equation that can be solved with standard methods, the solution requires some consideration. First of all, a shooting method cannot be applied, as one of the solutions of the homogeneous equation increases exponentially as X tends to infinity. Thus the two-point boundary-value problem has to be solved. Secondly, the possible appearance of very long, but very shallow, tails makes it advisable to choose the computational domain and the mesh size with care.
The size of the computational domain was chosen sufficiently large to comprise the expected exponential decay of the solution upstream and downstream, respectively, according to (32). The solutions are subject to strong upstream decay according to the root k 3 of (32). Since the wave crest is expected to be upstream of the bump in the cases to be considered, cf. Fig. 3b, the upstream boundary of the domain was chosen to be at X left = −20. In contrast, the downstream decay will be weak, as predicted by (36). Thus, the downstream extension of the domain ought to depend on the dissipation parameter β. In order to constrain the error, X right is determined from the condition β X right = 8, e.g. X right = 80 for β = 0.1 as a typical value of the dissipation parameter. The boundary conditions (29) are then implemented as H 1 = 0 at X left and X right , respectively. This approximation introduces errors of the order of 10 −6 in H 1 .
A second-order difference scheme for non-uniform grids was used to discretize the terms in (28). The grid point density near the ends of the computational domain was allowed to die out for the cases where a large domain is required. The derivatives H 1X and H 1X X X were discretized by a central 3-point stencil and a central 5-point stencil, respectively [3]. Equivalent one-sided second-order schemes were used at the domain boundaries. The finite difference formulas are given in Appendix E.
For the present cases, β values of the order of 0.1 are sufficiently large to allow the use of an equidistant distribution of grid points. The numerical solutions given below as examples were obtained with 12,000 uniformly distributed grid points. The computation was repeated with 6000 grid points in order to check whether the resolution was sufficient. The comparison of the computations with 12,000 and 6000 grid points, respectively, showed a difference of the order of 10 −6 in the height of the wave crest. The sparse nonlinear system of equations was solved with the Newton method. The Jacobian matrix was computed analytically. The solution for weak dissipation was used as an initial guess, if available, see Appendix F for details. The linear systems were solved with the LinearSolve command in Mathematica. The maximum value of the residuals after the Newton iterations was below 10 −8 in all cases.
In order to confirm that the solver gives correct results, a comparison with analytical results [47] is given in Fig. 4. Four cases of ramps are shown. The ramps have different non-dimensional heights, whereas the non-dimensional length and the dissipation parameter, respectively, are the same in all cases. The relative error of the uniformly valid asymptotic solution is of the order of β. The comparisons show reasonable agreement even for the rather large value β = 0.1. For considerably smaller values of β, the numerical and analytical solutions become indistinguishable from each other.  (29), in comparison to uniformly valid solutions [47]

Experiments
To perform experiments with the aim of verifying the theoretical results, a facility was available at Nihon University, Tokyo; see [23,53,54] for technical details. The channel width and length are 0.8 m and 15 m, respectively. A typical value of the volume flow rate per unit channel width isV ∼ 0.1 m 2 /s. With a channel width of 0.8 m and a water depth of approximately 0.1 m in the fully developed flow far upstream of the bump, the hydraulic diameter becomes d h ∼ 0.3 m, leading to Reynolds numbers Re, based on the hydraulic diameter, between 200,000 and 400,000. However, when the Reynolds number is determined for a particular experiment, it has to be taken into account that the viscosity of the water depends on the room temperature in the laboratory, which is subject to seasonal variations, among others.
In order to allow the application of bumps without ruining the facility, the channel bottom was covered with steel plates of 0.796 m width. Six steel plates were plane. Each had a thickness of 5 mm and a length of 1 m. An additional steel plate, containing a bump with constant cross section in lateral direction, was placed behind the third plate. The additional plate had the length of the bump; cf. Fig. 5. The shape of the bump cross section was an isosceles triangle. Observations for steep-slope channels [53] as well as estimates based on theoretical investigations on non-developed in-flow [33,35] indicate that placing three plane steel plates upstream of the plate containing the bump is sufficient to guarantee fully developed flow upstream of the wave, while the remaining three plane plates were placed downstream of the bump plate in order to simulate a plane bottom of large extent. The first steel plate was located 5 m downstream of the upstream end of the channel. The flow dropped over a head tank with 1.2 m height at the upstream end of the open channel. As the channel was sloped until the required steep slope, the supercritical flow was formed without using a sluice gate. The flow was a developing flow upstream of the location of the first steel plate. Along the three steel plates upstream of the bump, the flow underwent a transition from gradually varying flow to uniform flow with a depth of about 9.08 cm for Fr = 1.08. By measuring the water depth it was confirmed that the flow depth was constant upstream of the bump, indicating that the flow was indeed fully developed.
Results for two bumps will be presented below. One bump had a half-length of 195 mm and a height of 1.5 mm, the respective values of the second bump were 135 mm and 3.0 mm. For technical reasons it was not possible to use bumps with heights considerably smaller than 1.5 mm. The larger height of 3.0 mm was then chosen for the purpose of observing possible effects of the bump height. That resulted in non-dimensional bump heights ψ(0) of 23.3 and 46.7, respectively. Unfortunately, those values are rather large, whereas they are assumed to be of the order of 1 in the asymptotic analysis.
Concerning ramps, it was not possible at the time being to perform experiments to a sufficient extent and with sufficient accuracy. Future work may be devoted to that problem.
Based on the values given in [59], Table 17.1, the roughness of seamless steel plates can be estimated to be between 0.05 and 0.1 mm. To smooth the surface, it was covered with an oil film. Thus, the surface can be considered as hydraulically smooth. Knowledge of the hydraulic roughness makes the planning of the experiments easier, but it is of no relevance for the comparison with the analysis, provided the flow is fully developed far upstream of the stationary wave. This is accomplished by arranging the steel plates as described above, and by adjusting the channel slope α according to the upstream flow parameters. Then, the channel slope becomes of the order of α ∼ 3 × 10 −3 .
The water surface elevation was measured with a point gauge with ± 0.1 mm reading. The uncertainty of the measurement is about ± 1 mm. Apart from the turbulent fluctuations, the water surface was found to be practically steady, i.e. in accord with the assumption of steady mean flow. Video 1, which is available online, lends support to that observation. The video shows the water flow over a bump of isosceles triangular cross section with half-length 135 mm and height 3.0 mm; the bottom slope is 1/327, while the Froude number and the Reynolds number (based on hydraulic diameter) are Fr = 1.08 and Re = 400,000, respectively.
In order to check whether the surface elevation is also in accord with the assumption of two-dimensional flow, the surface elevation at different lateral positions was measured; see Fig. 6 for a typical result. It appears that the secondary flow, which is certainly present [57], has little effect on the surface elevation. Nevertheless, in order to minimize possible sidewall and secondary-flow effects on the measurements, laterally averaged values are used for the comparison with theoretical results.
Another assumption requiring justification concerns the viscous wall layer at the bottom. According to Sect. 3.2, first paragraph, the thickness of the viscous wall layer is assumed to be much smaller than the height of the bump. Using (4) for the friction velocity, the thickness of the viscous wall layer is estimated to be of the order of ν/ gαh r . In the experiments, the room temperature in the laboratory varied between 10 and 27.3 • C. Accordingly, the kinematic viscosity of the water varied between 1.31 × 10 −6 and 0.754 × 10 −6 m 2 /s. This gives a thickness of the viscous wall layer of about 24 µm and 14 µm, respectively. In view of those values, the bump heights of 1.5 mm and 3 mm, respectively, are well above the required minimum.
Comparisons of measured surface elevations with the theoretical predictions according to the present analysis are given below. The values of the dissipation parameter β were determined from (19) with the bottom slope α as measured at the facility for the particular experiment. The results to be shown have been selected from several experiments with the aim of verifying the existence of two stationary wave types. The experiments were performed in the well-established way [23,53,54], i.e. by adjusting the channel slope and the discharge in order to settle the required Froude number and the required supercritical flow depth. When the water surface profile was recorded transversely and longitudinally under steady flow conditions, which took about 5 min, a surprising observation was made. In contrast to what has been observed for a plane bottom with a region of enlarged roughness [62], the stationary surface elevation over the bump did not come close to the predicted solitary wave. The stationary single wave was rather of the type of the solutions of the second kind; see Fig. 7 for a typical example. Figure 7 shows that the position and the maximum height of the observed stationary wave are in good agreement with the predictions, whereas the width of the observed wave is somewhat larger, and the tail is a bit less pronounced, than according to the numerical solution of the second kind. However, the deviations are within the approximate error bounds, cf. the estimates given below. Far downstream of the stationary wave, after the flow had already reached a region of nearly undisturbed water depth, another stationary surface elevation was observed, as the measured value near X = 10 indicates. The latter surface elevation was apparently caused by a very small step in the bottom near the sidewall. Regrettably, it was not possible to remove that unevenness of the bottom in course of the experiments, but, as the flow is supercritical, an effect on the single wave quite far upstream of the unevenness can be ruled out.
In course of attempts to realize the stationary wave of the first kind (solitary-wave type), the following start-up process was discovered. A point gauge was inserted into the water down to the bottom and kept there for some time, e.g. 40-60 s. As a consequence, the wave resembling the solution of the second kind moved upstream of the bump, while the amplitude of the wave became larger. After the movement of the wave, the wave formation was stable when the point gauge was removed from the water. Figure 8 shows the results of measurements for a stationary wave of the solitary-wave type that was obtained in that way.
It is thought that those observations clearly indicate the existence of two types of stationary single waves, as theoretically predicted. Note that tails are visible in all results of measurements as given in Figs. 6, 7 and 8. With regard to quantitative comparisons of experimental data and theoretical results, the following error estimates can be made. The derivation of the extended KdV equation is based on an asymptotic expansion for Froude numbers close to 1 and large Reynolds numbers. Thus, the relative errors of the extended KdV equation are of the order of ε and 1/ln Re τ , respectively. With ε ≈ 0.05 according to (17) and 1/ln Re τ ≈ 0.12 according to (6), the error of the solutions of the extended KdV equation could be expected to be of the order of 10%. However, a particular problem is the bump with the rather large height of 3.0 mm. According to the asymptotic expansion described in Sect. 4, the forcing term in the extended KdV equation (28) gives rise to a relative error of the order of εβψ(0). With the values given in the caption to Fig. 8 one obtains εβψ(0) ≈ 0.20, i.e. an error of about 20%. Regarding the experimental data, the error in the water surface elevation measured with the point gauge is about h ∼ 1 mm, as already mentioned above. Using Eqs. (5a), (22) and (17), the error can be expressed in terms of the perturbation of the surface elevation as H 1 ∼ 3 h /2(Fr − 1)h r . With the parameters given in the captions to the Figs. 7 and 8, one obtains H 1 ∼ 0.2. Accounting for the error estimates for both theory and experiment, the deviations of the measured values of the surface elevation from the predictions appear reasonable. With regard to the numerical solutions it may be worth mentioning that it was not possible to find a numerical solution for the stationary wave of the second kind in case of the larger bump height, i.e. with the parameters as given in the caption to Fig. 8. Besides, in this case the numerical solution of the solitary-wave type exhibits a hump and differs considerably from the first-order analytical solution (35), i.e. the classical solitary wave, which is shown for the purpose of comparison in Fig. 8. Those peculiarities seem to indicate limits of applicability of the asymptotic expansions that form the basis of the present analysis, in particular with respect to the height of bumps.
If the predictions shown in Fig. 8 had been known a priori, a better spatial resolution of the experimental data could have been obtained locally at, and near, the predicted hump or other particular features of the numerical solutions. However, the numerical results were not available prior to conducting the experiments. As the results are very sensitive to the parameters, i.e. (Fr − 1), bottom slope and reference surface height, it was not possible to adjust all parameters in a particular experiment to a parameter set prescribed according to a numerical solution. Instead, the experiments were performed with parameters as close to the desired values as possible, and afterwards the numerical solutions were obtained with the data sets of the experiments for the purpose of comparison.
Finally, it might be worth mentioning that flow separation was not observed in any of the experiments with bumps. Obviously, the inclination angle of the bump flanks is too small to give rise to noticeable separation.

Conclusions
The present paper concerns turbulent open-channel flow over bumps or ramps of very small height. In the main theoretical part, a double asymptotic expansion was performed. Firstly, an asymptotic expansion for very large Reynolds numbers or, equivalently, very small slopes α of the channel bottom, could take advantage of splitting the flow field into two regions, i.e. the universal viscous wall layer at the bottom and the defect layer. For the latter, matching conditions with respect to the wall layer and boundary conditions at the free surface had to be satisfied. The novel feature in that approach is the treatment of the boundary conditions at the bumps or ramps. Secondly, in order to describe near-critical free-surface flow, the Froude number was written as Fr = 1 + (3/2)ε, and an asymptotic expansion was performed with ε as a small parameter. Though the governing equations of the defect layer and the dynamic boundary conditions at the free surface contain Reynolds stresses, turbulence modelling could be avoided by a proper choice of the relative orders of magnitude of α, ε and the non-dimensional height of the bump or ramp, B.
The results indicate that a very small unevenness in the channel bottom produces much larger surface elevations in steady near-critical turbulent open-channel flow. Shape and position of the surface elevation are governed by the steady-state version of an extended KdV equation. The extension consists of two parts: one is a linear term due to dissipation, which is known from previous investigations; the other one is a new forcing term that describes the effect of the unevenness in the bottom. For a ramp of a particular shape, the shape of the stationary surface elevation is found to be that of the classical solitary wave of inviscid flow over a plane bottom. The same solution is obtained in the limit of vanishing dissipation parameters for ramps and bumps of any shape. However, depending on the shape and size of the ramp or the bump, the position of the stationary solitary wave varies. With the exception of a ramp with a non-dimensional height of a particular value (eigenvalue), a long, shallow tail is attached to the solitary wave. For finite values of the dissipation parameter, stationary single waves resembling solitary waves with tails were obtained by numerical solutions of the extended KdV equation (28) with boundary conditions (29).
In addition to the waves of the type of solitary waves, stationary single waves of a second kind were obtained as solutions of the extended KdV equation. Their amplitudes are considerably smaller than those of the solitary-wave type. All solutions of the second kind display a long tail.
The asymptotic analysis gives results that are exact in the limit of infinite Reynolds number Re, Froude number Fr approaching the critical value 1 and non-dimensional bump or ramp height B approaching zero. In experiments, however, one deals with finite values of Re (or α), ε = (2/3)(Fr − 1), and B. That limits the accuracy of the theoretical results for the purpose of comparison with measurements. The accuracy of the measurements, on the other hand, is also rather limited, as the water surface elevation is relatively small. However, there can be no doubt that the experiments prove the existence of two fundamentally different steady states. In particular, it appears to be the first time that a solution of the second kind of an extended KdV equation has been observed in an experiment. Position and maximum surface elevation of the stationary wave of the second kind were predicted with good accuracy, without applying empirical constants. For the stationary waves of the first kind, i.e. those resembling stationary solitary waves, the agreement between experimental data and theoretical predictions is less good, but still within the estimated error bounds.
Two steady states are not uncommon in fluid mechanics. Examples that are close to the present one have already been mentioned in the introduction in the discussion of inviscid-flow solutions. Often, different start-up processes lead to different steady states. In the present case, it is the insertion of the point gauge that causes a start-up process that differs from that without an initial disturbance. In Appendix D, single waves of the solitary-wave type, varying slowly with time, are considered for the case of flow over a bump. That allows interesting, though limited, conclusions regarding stability of the stationary single waves of the first kind. Investigations on the stability of the stationary single waves of the second kind were beyond the scope of the present work. Clearly, the stability of the stationary single waves described in the present paper deserves further investigations.
In the theoretical part of the paper, equations were derived that relate the most important flow quantities, such as velocity components, pressure and shear stress at the bottom, to the free-surface elevation, which is to be determined as a solution of the steady-state version of the extended KdV equation. Measuring other flow quantities in addition to the surface elevation, though desirable, was not possible with the facility available for the present experiments. It may be noted, however, that the Reynolds shear stress distribution at the location of the first crest of an undular jump was measured by Lennon and Hill [41]. A comparison with solutions of (28) and ψ X ≡ 0 showed reasonable agreement [35]. In view of the similarity between the first wave of an undular jump and the solitary wave, this may be considered as lending further support to the present analysis.
Concerning future work, experiments with ramps are desirable. Preliminary investigations led to promising results, but they are still incomplete. Thus it was decided to refrain from including them in the present paper. It will also be of interest to compare the analytical results with numerical solutions of the full equations of motion with modelled Reynolds stresses, e.g. by applying the iteration method described in [61]. Work in that direction is in progress.

A Justification of using the logarithmic law of the wall
According to the generalized law of the wall [59], p. 547, a non-dimensional coupling parameter K characterizes the effect of the pressure gradient on the law of the wall. For K → 0, the generalized law of the wall reduces to the logarithmic law. In the present notation, K is defined as Introducing non-dimensional variables according to (2), (5a) and (6) gives Substituting for δ and expandingP according to (1) and (25), respectively, and introducing β with (19), one For the parameters of the experiment, K-values of the order of 10 −3 are obtained, justifying the application of the logarithmic law.

B Derivation of the extended KdV equation (28)
The following derivation of the extended KdV equation (28) follows [60], accounting for different boundary conditions at the bottom due to the presence of bumps or ramps. The dependent variables are expanded in terms of powers of ε as shown in (22) for the non-dimensional surface elevation as an example. To keep the analysis free of turbulence modelling, several important points have to be observed. First, the volumetric mean velocity of the fully developed flow far upstream is taken as the reference velocity, as already noted in Sect. 2. Secondly, the surface height of the fully developed flow, i.e. the reference heighth r , ought to be known for given values of bottom slope, α, and volume flow rate,V . This implies that the friction coefficient of the channel, c f r , is known, as the force balance gives h r = (c f rV 2 /2gα) 1/3 . Thirdly, the non-dimensional "velocity defect" U = U (Y ) is introduced such that the local value ofŪ in the fully developed flow far upstream differs from its volumetric mean, i.e. 1, by the amount √ α U ; cf. the comments below (28). In addition, the perturbation of the defect velocity profile, which leads to a perturbation ofC as defined by (13), will appear only in terms of higher order than presently considered. For that, it is of importance to apply the logarithmic "law of the wall" ("overlap law") in a suitable form, such as (11). On that basis, the asymptotic analysis can be performed as follows.
With the basic state given by the simple relations the expansion of the continuity equation (7) leads to U 1X + V 1Y = 0, which can be integrated to obtain where the expanded boundary condition at the bottom, Eq. (10), has already been satisfied. Before expanding the momentum equations (8a) and (8b) for small values of ε, the small parameters δ and α are expressed in terms of ε according to (1) and (18), respectively. Expanding then (8b) gives P 1Y ≡ 0, i.e.
Then expanding (8a), one obtains U 1X = −P 1X . This can also be integrated. The free function of integration is determined from relations characterizing the fully developed flow far upstream, i.e. P 1 ≡ 0 , is the non-dimensional velocity defect, which is introduced as described above. Therewith one obtains the relation Next, the following relations follow from the dynamic boundary conditions (16b) and (16a), in this order: Taking (52) and (53) into account, the integration in (50) can be performed, with the result Finally, the logarithmic law (11) is expanded together with (13). Concerning the expansion of the latter equation it suffices to formally writeC(X ) =C r + εC 1 (X ), whereC r denotes the value ofC in the reference state, i.e. for the fully developed flow. Using (6) and (12), one obtains the following relation: (56) Far upstream, the velocity perturbation is given by the velocity defect of the fully developed flow, and the perturbation of the Reynolds shear stress vanishes, i.e.
Comparing (57) with (56) shows that the logarithmic terms in the latter equation may be substituted according to the relation √ α to obtain where Fr has been replaced by 1 in leading order, and √ A U (1) has been eliminated with the help of (52) and (53). Note that Y = 0 refers to the channel bottom in terms of defect-layer coordinates, i.e. −(U V ) 1 (X, 0) is equal to the first-order perturbation of the bottom shear stress.
With (52)- (55) and (59), the results of the expansion up to first order are complete. Obviously, either H 1 (X ) or P 1 (X ) remain free in the framework of the first-order equations. Thus, the second-order equations have to be inspected for solvability.
Expanding, first, the momentum equation (8b) up to second order, making use of the first-order results given above, integrating with respect to Y and determining a free function of integration from the dynamic boundary condition (16b) at second order, gives where the last term on the right-hand side of the equation is a function of Y only. Secondly, the momentum equation (8a) is also expanded up to second order, retaining the term of the order βε 2 , with β defined in (19), to obtain U 2X may be replaced by −V 2Y according to the expanded version of the continuity equation (7), and, furthermore, P 2 may be eliminated using (60). The equation hence obtained can then be integrated with respect to Y, accounting for the boundary condition at the bottom, Eq. (10). Introducing, once more, the first-order results and observing that, by definition, one obtains the following relation: On the other hand, the second-order kinematic boundary condition, which follows from (15) upon expanding, gives i.e. (28). Note that U (1) has cancelled, so that knowledge of the velocity defect is not required.

C Conservation of momentum flow
As the steady flow far upstream (X → −∞), as well as far downstream (X → +∞), is assumed to be fully developed, the momentum flow rate far upstream is the same as far downstream, and the forces acting on the total volume of fluid must balance. The hydrostatic pressures far upstream and far downstream balance each other. Thus, the forces that remain to be taken into account in the over-all force balance are the X -component of the gravity force acting on the liquid, the force due to the pressure acting on the bump or ramp, and the bottom friction force. Since fully developed flow is characterized by an equilibrium of forces, only the forces in excess to those acting in the fully developed flow have to be considered. This gives the following force balance: The coefficients in front of the first and the third integral cancel according to (4). Substituting for δ, B and U 2 τ according to (1), (20) and (12), respectively, and introducing the expansions (22) to (26) gives, to the first order, With the integral boundaries (21a,b) for ψ, one finally obtains the integral relation (30).

D Slowly varying solitary wave in turbulent free-surface flow over a bump
For turbulent free-surface flows with ensemble-averaged quantities depending on time t, the analysis of Sects. 3 and 4 can easily be extended to obtain the following extended KdV equation: where T is the non-dimensional time variable, which is defined by the equation Applying a method due to [63], pp. 294-296, to the bump with isosceles triangular cross section defined in Sect. 5.2.3, one obtains the following expression for the slowly changing surface elevation: with the evolution equation for the slowly varying wave speed V (T ). When evaluating the integrals that lead to the right-hand side of (71), The stencil used at the right domain boundary point i = n is left-biased accordingly.

E.2 Third derivative
At the left boundary point i = 1, a right-biased 5-point stencil is used For the point next to the left domain boundary i = 2 For the interior points i = [3, n − 2], a central 5-point stencil is used, cf. [3] Accordingly, left-biased stencils are used for i = n − 1 and i = n.

F Initial guess for the Newton iteration
In the numerical solver described in Sect. 6, the following initial guesses are used to start the Newton loop. For the numerical solutions in Fig. 4, the uniformly valid solutions from [47], as shown in the figure, are used as initial guesses. For the bump in Fig. 7, a uniformly valid analytical solution of the first kind is not available yet. Instead, we selected  (1). For bumps, ψ ∞ = 0. We chose X (0) m = −2 and X (0) m = −1 for the stable and unstable solitary wave, respectively. For the solution of the second kind, the first-order inner solution according to [49], Eq. (6), was used as the initial guess. For the result in Fig. 8, the displayed classical solitary-wave solution was used.