Hydrodynamic forces on emergent cylinders in non-uniform flow

The main objective of conducting numerical simulations of flows in rivers with vegetation is to investigate the complex flow dynamics involved in non-equilibrium conditions. In such cases, it is inappropriate to apply the drag coefficient CD, which is typically derived based on uniform flows involving groups of infinitely long cylinders. This paper presents a method for evaluating the drag forces acting on emergent obstacles for non-uniform open-channel flows. This method is devised based on two sets of experiments: on flows with small-diameter cylinders, focusing on the water surface profiles through the group; and on flows with large-diameter cylinders, focusing on the local pressure distribution and local water surface profile around a target cylinder. In addition to the conventional drag force expression that includes CD, two new terms are proposed to account for the effects of water surface variation and pressure gradient in non-uniform open-channel flow on the drag. The first of these terms, which introduces the use of the Froude number to account for the effect of water surface variation, is derived theoretically and evaluated against past and present experimental results under uniform-flow conditions. On the other hand, the second of these terms, which includes the representative length of the separation zone to evaluate the effect of pressure gradient, is confirmed to be a necessity through numerical calculation of the longitudinal water surface profile in emergent cylinders. The incorporation of these two terms using a simple unified expression can help improve the accuracy of numerical simulations for practical problems of flows with emergent obstacles. Conventional drag-force terms are shown to be insufficient for calculations of water flows with emergent obstacles. The effect of water surface variation on the drag force is derived theoretically as a function of the Froude number. For accelerating flows, additional drag forces are revealed by local pressure distribution and water surface profiles. Conventional drag-force terms are shown to be insufficient for calculations of water flows with emergent obstacles. The effect of water surface variation on the drag force is derived theoretically as a function of the Froude number. For accelerating flows, additional drag forces are revealed by local pressure distribution and water surface profiles.

(1) Conventional drag-force terms are shown to be insufficient for calculations of water flows with emergent obstacles. (2) The effect of water surface variation on the drag force is derived theoretically as a function of the Froude number.

Introduction
The effects of vegetation on water flow, turbulence, and sediment transport have attracted considerable interest from researchers in a variety of fields [1][2][3][4]. For example, it has been demonstrated that the characteristics of canopy flows should be clearly determined to properly understand the behaviors of wind flows involving forest canopies [5,6] and groups of buildings in urban areas [7], as well as the behaviors of water flows involving coral reefs [8] and in gravel-bed rivers with boulders [9][10][11]. It is also important to know the resistance characteristics of a group of cylinders in fluid flow because of the occurrence of this structure in a number of applications of fluid dynamics, such as heat exchangers, where arrays of cylinders are a common configuration [12]. Generally, in flow calculations, multiscale obstacles are divided based on two characteristic scales: large-scale topography, which can be resolved using a computational mesh explicitly as a boundary shape; and small-scale topography, which is evaluated using subscale models [9,11]. In a number of studies, numerical simulations of flows through vegetation have been conducted using fine meshes to resolve the effects of obstacles [7,13]; in others, the effects of vegetation have been modeled using temporal-spatial averaging of the governing equations [3,5,6,8,[14][15][16][17][18]. Based on the complexity of vegetation conditions in rivers, the second aforementioned approach is considered to be suitable for practical applications [18][19][20]. When averaged over space in Reynolds-averaged Navier-Stokes equations for local flows, in addition to the dispersion term due to sub-average-scale flows, the fluid force terms appear in conjunction with Gauss' divergence theorem in the momentum equation integrating the local pressure distribution and shear stress acting on the objects [6,15]. The problem, as with the closure problem for turbulence, is how to represent the hydrodynamic force terms induced by sub-average-scale flows around the boundaries and distinguish them from the average-scale flow [21]. The drag-force terms in the momentum equations are evaluated using drag coefficients, which determine the accuracy of the models [22]. The drag coefficients also appear in equations for evaluating the turbulence intensity and dispersion coefficients of passive particles in water flows [3,14,17,23,24].
Meanwhile, several studies have contributed to determining the drag induced by vegetation [25][26][27][28][29]. However, although actual vegetation tends to be irregularly shaped and flexible [18], circular cylinders have been employed to model rigid vegetation in calculations regarding the characteristics of canopy flow. The drag coefficient C D of a circular cylinder in a free-stream flow (defined as a uniform flow that lacks a free surface, and typically distinguished from a uniform open-channel flow) is known to be a function of the Reynolds number R = Vd/v (V, free-stream velocity; d, cylinder diameter; v, kinematic viscosity); the same is true for a sphere among several representative immersed bodies [30][31][32]. The variation in the drag coefficient of an array of cylinders depends on the space between the cylinders and the density of the array. Several researchers have contributed to correlating drag coefficients of cylinder arrays with well-established relationships between the drag-force coefficient and Reynolds number for a single cylinder. For example, Koch and Ladd [25] and Tanino and Nepf [26] applied the equation proposed by Ergun [33], which was devised for pressure drops in packed columns, to the drag coefficient and Reynolds number function. Kothyari et al. [34] proposed a relationship between the drag coefficient and Reynolds number for subcritical flows as a positive logarithmic function of the density, showing that the drag coefficient increases rapidly for low density values. Cheng [35] proposed a formula for the drag coefficient of a cylinder array by expanding the relationship for an isolated cylinder based on the concept of a pseudo-fluid model [36] and validated the formula against several drag-force coefficient datasets.
On the other hand, in contrast to most canopy flows, open-channel flows, such as in rivers, involve the Froude number F = V/(gh) 0.5 (g, gravity; h, water depth), an important dimensionless parameter that influences the drag coefficient, especially for emergent obstacles in shallow-water flows. In fact, the drag coefficients of hydraulic structures are related to the Froude number [31]. However, only a few studies have investigated the effect of the Froude number on the drag coefficients of emergent cylinders in open-channel flows. Furthermore, Kothyari et al. [34] indicated that, whereas the drag coefficient is not substantially influenced by the Froude number in subcritical flows, the coefficient decreases with respect to the Froude number in supercritical flows. The results of a study by Huai et al. [37], in which a large-eddy numerical simulation model was used in conjunction with three identical-density conditions, showed that C D largely decreases as F increases. It is significant that C D is smaller in supercritical flows than in subcritical flows; however, the underlying mechanism for this distinction has not yet been elucidated.
Most previous studies on the evaluation of drag forces and flow resistance were conducted based on uniform open-channel flows. However, to fulfill the objective of applying drag-force coefficients to numerical simulations of flows in natural rivers [19,20], a method for evaluating the effects of non-equilibrium flow in average scales on the resistance must be developed. For example, Busari and Li [38] investigated the interference effects of finite-length vegetation in gradually varied flows, whereas Uchida et al. [11] proposed a dynamic rough-wall law to account for the effect of non-equilibrium flow near a bed with submerged boulders on the flow resistance. However, despite these efforts, there has been insufficient research on methods for evaluating the direct effects of non-equilibrium open-channel flows on the drag force.
To take into account the resistance of vegetation in rivers, it is necessary to elucidate the effects of variations in the water surfaces around obstacles. For non-equilibrium flows with pressure gradients in an average scale, because the separation regions behind the obstacles exhibit low flow velocities [12,39,40], pressures acting on the backs of the obstacles are inferred to be affected by the pressure downstream in the scale of the separation length.
The objective of this study is to develop a method for evaluating the drag coefficient of emergent cylinders in an array in a non-uniform open-channel flow. The proposed method is to be designed for numerical simulations. In addition to the conventional drag-force expression, which includes the drag coefficient, this paper proposes two new terms for evaluating the effects of water surface variation in a local scale and pressure gradient in an average scale on the drag force for emergent cylinders. For the first term, a drag-coefficient equation is derived theoretically as a function of the Froude number. This equation is compared with drag coefficients from previous experiments and our laboratory experiment on uniform flows without the introduction of any calibration parameters. For the second term, the representative length of the separation zone is introduced to evaluate the effect of the pressure gradient in the average scale on the drag force in non-uniform open channel flows. The calculation of gradually varied flows with emergent cylinders reveals the necessity of the term to account for the effect of the pressure gradient on the drag forces in accelerating flows. The validity of the term is supported by the results of supplementary experiments that focus on the distribution of pressure on a cylinder in a group. Based on this pressure distribution, the negative pressure gradient for accelerating flow decreases the pressure in the separation zone behind the cylinder. It is then demonstrated that several non-uniform water-surface profiles, including accelerating and decelerating flows, can be reproduced well using a simple unified expression incorporating these two terms.

Experimental setup
Two sets of experiments were conducted in this study. The first set of experiments investigated longitudinal water surface profiles using an experimental channel, in which several small-scale emergent cylinders were installed, and focused on developing a method for evaluating the resistance. To clarify the effects of water surface variation and pressure gradients on drag force, three types of flow behavior, i.e., uniform, accelerating, and decelerating flows, were examined. The second set of experiments involved a supplemental study that elucidated the mechanism of the hydrodynamic force and pressure distribution on emergent cylinders. This experiment was conducted to measure the local water surface variations around large-scale emergent cylinders and the distribution of pressure acting on a cylinder in the group.
The first set of experiments was conducted in a water-recirculating flume (length, 3.0 m; width B, 0.3 m) with varying channel slopes. The flows were adjusted to either uniform or non-uniform flow using a weir installed downstream. For this study, we regard uniform flow as that when the water surface gradients coincide with the channel gradient. The emergent cylinders were circular wooden dowels (diameter d, 10 mm; height, 20 cm). The Styrofoam floor of the flume was drilled, and the dowels were installed individually in longitudinal and lateral intervals Δx and Δy, respectively, of 0.0344 m in a staggered array from the downstream end of the channel ( maximize the flow resistance. The density of cylinders per unit area, λ, is defined as λ = Nπd 2 /4, where N is the number of cylinders per unit area (present experimental conditions: N = 1/(ΔxΔy)), and d is the cylinder diameter. The density λ = 0.067 in this experiment is in the range of actual vegetation densities in the flood plains of Japan [41].
The experiments consisted of Case L and S series: In the Case L series (Fig. 1b), the cylinders cover the last 2.0 m of the flume, whereas in the Case S series (Fig. 1c), the cylinders cover only 1.0 m of the flume between the upstream and downstream ends. Flow depth was measured using a point gauge with a 0.1-mm vernier scale mounted on a carriage that could be moved along the flume rails. The water depth h was measured along the longitudinal axis on the centerline of the flume at intervals of 0.2 m. The resistance evaluation in this study is based on the results of this water surface profile measurement, i.e., the focus of this study is the modeling of phenomena that can affect the larger-scale analyses of one-dimensional measurement approaches, such as water surface profiles. On the other hand, the effects of the upstream and downstream ends of a cylinder group on the internal structure, such as on the velocity distribution and turbulence in two-or three-dimensional flow structures, are not included in this study. Investigations regarding these aspects of flow behavior will require further advancement of the results of this research. The volumetric discharge was measured using a triangular weir with a tank 0.8 m in length, 0.4 m in width, and 0.4 m in height. The average cross-sectional velocity (bulk velocity) U was defined using the mean depth h as U = Q/(hB), where Q and B are the discharge and width of the flume, respectively.
There are several possible definitions of representative velocity for evaluating canopy drag, including bulk velocity U [16-18, 42, 43], pore velocity U p = U/(1 − λ) [18,26,34,41], and constructed velocity U c = U/(1 − d/Δy) [44]. The undisturbed velocity, which is defined as a local velocity removing the target obstacle, is considered appropriate to evaluate a local pressure change by installing the target obstacle [45][46][47]. Because the local pressure change around the obstacle induces the hydrodynamic force, the representative velocities include the bulk, pore, and constructed velocities, which can be regarded as attempts to estimate the undisturbed velocity. However, the undisturbed velocity is affected by other cylinders and hydraulic conditions, and thus a proper method for the evaluation of the undisturbed velocity has not yet been established [47]. The definition of representative velocity is also relevant to discussions regarding the representative flow velocity for the Froude number, which is relevant to deriving the depth profile equation (see Eq. 8). In this study, the vegetation density was relatively low, and its effect was small [18]. Furthermore, as in many numerical calculations on flows involving vegetation [16,17,20], the bulk velocity was employed to focus on the effects of the water surface and pressure gradient on the drag coefficient. Nevertheless, there remain questions regarding the representative velocity for future studies to resolve. Table 1 summarizes the bed slope and discharge conditions for the experiments on uniform, accelerating, and decelerating open-channel flows. For Cases L1-L6, the differences in the gradients between the water surface and bed, calculated using measurements obtained at 0.2-m intervals over 1.4 m, were less than 0.4 × 10 −4 . Large values of R were used to decrease the dependence of the drag coefficient on the Reynolds number [26-29, 34, 35, 48]. To distinguish the effects of the Froude number from those of the Reynolds number, a more detailed investigation will be necessary because the effect of the Reynolds number on the drag coefficient varies with cylinder conditions [25][26][27][28][29]. Thus, instead of analyzing a large amount of drag coefficient data, we examine the effects of the Froude number and water surface theoretically in the following section.
The second set of experiments, a supplementary study using larger-scale cylinders (diameter 0.  Table 2. In this set of experiments, the representative velocity for calculating the pressure coefficient and drag coefficient was measured using a MicroADV (SonTek) probe based on an undisturbed flow, i.e., without the cylinder of interest, to eliminate the effect of high velocities along the sidewalls, which were due to the small number of cylinders in the transverse direction (each row comprised only two or three cylinders, as shown in Fig. 2). The vertical distribution of pressures acting on the center cylinder of the group was also obtained. Figure 3 shows the pressure-measurement cylinder installed in the center of the group (x = 0 in Fig. 2). This 10-cm-diameter cylinder, along the side of which six equidistant holes were drilled (3 cm apart), was used to measure the pressure at different water depths. The cylinder was equipped with 5-mm-diameter tubes, which were connected to manometer glass tubes with a diameter of 0.6 mm. The pressure p on the cylinder was measured at an arbitrary angle θ with respect to the flow using the manometers as the pressure-measurement cylinder was rotated. The measurements were performed at intervals of π/18. The accuracy of the pressure measurement using the manometers was estimated to be 0.2-0.3 mm based on averages over the vertical pressure points for each angle.

Equation of water surface profile for flow with cylinders
To determine the relationship between the water surface profiles and drag forces on the vegetation, we first consider the momentum equation in the streamwise direction for shallow-water flows (hydrostatic pressure and no vertical velocity distribution) in a wide rectangular channel: where f t is the external force acting on the water column in a unit area, including the shear stress on the channel bed and drag forces on the vegetation. For a steady-flow condition, we have where q is the unit flow discharge. By substituting Eq. (2) into Eq. (1), we obtain where F 2 = q 2 /gh 3 . For the case involving no cylinders, f t = τ 0 = ρgS f (τ 0 , bed shear stress; S f , friction slope); therefore, the well-known governing equation for gradually varied flows [30,49] can be obtained. For flows with emergent-cylinder vegetation models, the total force f t includes a term for the drag force acting on the vegetation: where Nf d is the drag force acting on the cylindrical obstacles in a control volume with a unit area, and f d is the drag force on a cylindrical obstacle. Substituting the Manning roughness coefficient n in SI units and the drag coefficient C D , we obtain For evaluations involving bed shear stress, we did not consider deformations of the vertical velocity profile, as in many two-dimensional analytical studies on rivers [16][17][18][19][20]. It is known that the vertical velocity profile varies between accelerating and decelerating flows [50] and with the resistance by the cylinder group [44]. In addition, because the presence of the cylinder group results in horseshoe vortices, a secondary resistance due to the cylinder Pressure-measurement cylinder with manometers group, apart from the drag forces acting on them, is also generated [11]. To accurately evaluate the total resistance, an advanced model capable of calculating the non-equilibrium velocity profile [11] will be necessary. However, as indicated afterward, because the contribution of bottom shear stress to the total resistance is assumed to be much smaller than the drag force on the cylinders, the first term of Eq. (5) was simplified based on an assumption of a uniform velocity profile. When Eq. (5) is substituted into Eq. (3) for non-uniform flow, the water depth profile for a gradually varying flow with emergent vegetation can be obtained [38,49]: The results for uniform flow (dh/dx = 0) can then be deduced: For a flow that involves dense emergent vegetation, the drag coefficient term is known to be dominant over the bed friction term [26]. The bed friction term with n = 0.010 estimated for the channel bed material (Styrofoam) used in this study is much smaller than the drag force terms (< 1% for Cases L1-L6 outlined in Table 1). In this study, the drag force acting on the emergent cylinders in the uniform open-channel flows (Cases L1-L6) was evaluated using Eq. (7), using averaged gradients of the total head U 2 /2 g + h + z, instead of S 0 , to minimize the effects of non-equilibrium flow.

Drag in uniform open-channel flows
In this study, the drag f d , defined in Eq. (4), on an emergent obstacle in open-channel flows comprises three components: the base component f 0 , which excludes the effects of water surface variation and pressure gradients or drag forces acting on an infinitely long cylinder for deep uniform open-channel flows ( F = 0); and two additional components f s and f p , associated with variations in the water surface and the pressure gradient attributed to the drag force, respectively.
This section discusses f s ; the other component f p is discussed further in the paper in the section on the calculation of the water surface profile.
The base component f 0 has already been analyzed in several past studies. Typically, the drag force is composed of pressure and friction components, the former being dominant in flows with high Reynolds numbers [30,32], in which the drag is attributed directly to the relative pressure Δp acting on an object in the flow. Herein, the relative pressure on a vertical cylinder is expressed as Δp = p − p 0 , where p is the actual pressure acting on the cylinder surface, and p 0 is the pressure of the "free stream" or "undisturbed flow" [45][46][47] at that point without the objective cylinder. The pressure drag is often referred to as "form drag" because of its strong dependence on the shape or form of the object [30]. In general, pressure is described as a dimensionless pressure coefficient C p based on the representative velocity U under free-stream or undisturbed-flow conditions as follows:  8), is applied such that the other additional term is not considered (i.e., f p = 0). The water surface variations Δh are induced by the pressure variation Δp around the cylinder under the water surface. Given that negative gauge pressures are untenable, the water level is reduced to the height of the zero-gauge pressure. To derive the effect of water surface variation around an emergent obstacle on the drag coefficient, we considered a typical pressure distribution around a cylindrical obstacle at an arbitrary water depth h (Fig. 6a). When the hydrostatic pressure in the undisturbed openchannel flow ρg(h − z) (where z is the height of the section) is subtracted from the actual pressure p, the result is Δp, distributed as shown in Fig. 6b. Based on an assumption of a constant Δp in the vertical direction and in the hydrostatic pressure distribution, Δh(θ) is defined as where C p (θ) is the pressure coefficient. However, because of the range of variation in the water surface (z > Δp min /ρg, Δp min : minimum value of the relative pressure Δp), the effect of the water surface needs to be considered in the relative pressure distribution. Figure 6c shows the pressure variation at z = 0 where Δp = p, whereas Fig. 6d shows an example between the highest water level (in front) and lowest water level. In this range, a zero-pressure value occurs above the water surface, and the distribution of the actual Δp values is different from that shown in Fig. 6b. This induces an additional force as an effect of the water surface variation in open-channel flow.
To formulate the problem, we focused on the water surface variation component of the drag f s in Eq. (8). Based on sections at heights below the undisturbed water depth h, an additional pressure force Δf is required to compensate for negative pressures at angles with negative Δh, at which ρg(h − z) + Δp is negative. The additional differential force due to the water surface variation component, df s , is defined as Conversely, based on sections above the undisturbed water depth h, the same additional pressure as that appearing in Eq. (14) was applied under the assumption of a hydrostatic pressure distribution. Therefore, Eq. (14) is considered to be a definition of the infinitesimal force required to induce the water surface variation component of drag The relative drag coefficient C D /C D0 is proportional to the square of the Froude number, with the proportionality factor composed of C D0 and β Cp . The coefficient β Cp is defined as the second-order moment of the pressure coefficient around the cylinder, expressed as

Experimental results for uniform flows
The drag coefficient for base component of the drag force, C D0 , and the coefficient β Cp , calculated based on the distribution of the measured pressures on the cylinder, are listed in Table 2 and are compared with values calculated based on free stream-flow conditions [30]. Note that the coefficients β Cp for accelerating flow that are outlined in Table 2 include the incremental drag due to the pressure gradient, which is described further in the paper, and are defined differently from the C D0 based on Eq. (8). Figure 7 shows the variations in C D0 and β Cp with respect to vegetation density. Whereas the drag coefficient C D0 increases as the vegetation density increases, the coefficient β Cp tends to decrease. This is because the pressure coefficient C p used to determine β Cp using Eq. (19) depends on the Reynolds number and the conditions of the cylinder array. Moreover, the symmetry between the variations in the drag coefficient C D0 and the variations in the coefficient β Cp reveals a high correlation relationship between them. However, the available data are limited, and further study will be necessary to accurately predict the effect of water surface changes on the drag coefficient. Nonetheless, because all values of β Cp calculated herein, as shown in Fig. 7, are negative, it can be inferred that the drag coefficient in an open-channel flow decreases as the Froude number increases, which is consistent with the findings of previous studies [34,37]. Figure 8a shows the drag coefficient C D as a function of the Froude number F for the uniform-flow conditions listed in Table 1. Figure 8a enables comparisons between the experimental values in different hydraulic conditions and theoretical values based on Eq. (18). The drag-force coefficient based on the absence of water surface effects, or F → 0 , i.e., C D0 , is necessary for Fig. 8. The value C D0 = 1.22 is obtained to minimize the root-mean-square error (RMSE) between the measured drag coefficient and the drag  Table 2 is a result of a three-cylinder row array with high velocity along the sidewall, whereas C D0 = 1.22 when a sufficient number of cylinders in a row are used in the experiment. However, the differences in C D0 can be explained by the blockage effect, because the C D0 values 1.22 and 1.79 become 0.76 and 0.81, respectively, after correction for the blockage effect, i.e., C * D /C D = (1 −D/Δy) 1.35 (C * D : corrected C D for blockage effect) [51]. For comparison, Fig. 8 also includes previous experimental results from Kothyari et al. [34] and Ishikawa et al. [43] on emergent cylinders in uniform flows, in which the values of C D0 were obtained in the same manner as described previously, i.e., using Eq. (20). In line with the variations in C D0 and β Cp shown in Fig. 7, the experimental data are compared with calculation results based on possible ranges of values of C D0 and β Cp . Most of the experimental data are explained well within the range of the calculation results, except for a number of scattered values with low R . This result demonstrated how the influence of the Froude number on the drag coefficient came to be overlooked. However, according to the relationship between C D /C D0 and F (Fig. 8a) for higher values of F , the drag coefficient has an apparent tendency to decrease as the Froude number is increased, even in subcritical-flow conditions. This decrease in C D /C D0 with respect to increasing F was also confirmed when we compared the drag coefficients for Cases L3-L7 (as outlined in Table 1), in which F varied within a small range of R . The dependence of C D on R and F is discussed as follows.
In Eq. (18), the drag coefficient C D is used in conjunction with C D0 , β Cp , and F . The values of C D0 and β Cp are defined by the pressure coefficient C p , as expressed in Eqs. (12) and (19), respectively. Based on an assumption of a hydrostatic pressure distribution, as described earlier, C p is supposed to be independent of F , whereas C D is a function of F and C p . On the other hand, C p around a cylinder can be expressed as a function of R based on the similarity principle with respect to the Reynolds number for the Navier-Stokes equations [32]. From this, C D can be considered a function of the two independent variables F and R, as C D (F, R ). The empirical Ergun formula has been applied to the calculation of drag coefficient on aligned cylinders [25,26], as follows: Equation (20) represents the relationship between the drag coefficient and Reynolds number of a single body immersed in a free-stream flow with a relatively high Reynolds number [26]. The coefficients α 0 and α 1 of the equation have been investigated using experimental data and are known to be a function of the array characteristics.
If the C D0 in Eq. (18) is obtained using Eq. (20), and a constant β Cp is assumed, the value of C D can be expressed as a function of R and F as follows: Calculating the total derivative of dC D of Eq. (22) produces Eq. (23): Referring to Tanino and Nepf [26] for the order of α 0 , we investigated the dependence of C D on F and R for C D0 = 1 and α 0 = 100 (Fig. 9). Even in subcritical-flow conditions (log 10 (F ) < 0), the dependence of C D on F in flows with large Reynolds numbers is not 2 > 5 because the calculation produces a negative value for C D . More specifically, for the constraint condition that satisfies the positive-water-depth condition h + Δh > 0, the limitation on Eq. (18) is derived from Eq. (13) using C pmin , which takes the minimum water depth and value of C p as F 2 < − 2/C pmin (C pmin < 0), based on an assumption of free-stream C p (θ). For the condition F 2 > − 2/C pmin , the bed behind each cylinder was dried. Equation (13) for the Δh behind the cylinder must be modified to maintain the water depth at h + Δh > 0; thus, a penalty function is applied: Equation (24) can be introduced to satisfy the constraint condition h + Δh > 0 for the calculation of β cp using Eq. (19), in conjunction with the measurement of the pressure coefficient distribution based on a uniform open-channel flow ( Table 2). At high Froude numbers (F 2 > − 2/C pmin ), the coefficient β Cp must increase when the C p ′ in Eq. (24) is used as the C p in Eq. (19). For this uniform open-channel flow, the value of β Cp increases from a constant value (β Cp = − 0.24, − 0.40, − 0.51) for the supercritical-flow condition, eventually becoming positive at high Froude numbers (Fig. 10), indicating that the drag force on the emergent object will increase at high Froude numbers. Based on the variation in the coefficient β Cp , the theoretical values for the drag coefficients are modified, as indicated by the dashed lines in Fig. 8b. The relative drag coefficient C D /C D0 attains its minimum value and gradually increases as the Froude number is increased. However, the theoretical lines based on Eq. (24) (dashed lines in Fig. 8) still underestimate the measured drag forces, which deviate from the solid theoretical lines over C D0F 2 = 1, where there is still little discrepancy  between the solid and dashed lines. This discrepancy between the experimental results and the modified theoretical values is attributed to the occurrence of a vertical distribution of pressure change Δp induced by the vertical velocity distribution. By contrast, a uniform velocity distribution and a hydrostatic pressure distribution were assumed in the derivation of Eqs. (18) and (24) even for high Froude numbers. Because the velocity near the water surface is higher than the depth-averaged velocity, the actual effect of the water surface on the drag force is considered to be greater than that assumed in the present theory. Thus, the assumption that C p is independent of F , which was used to introduce Eqs. (18) and (21), is no longer satisfied at high F conditions. However, although detailed investigations will be necessary for more accurate estimations of the drag coefficient, particularly for flows with high Froude numbers, Eq. (18) can be adopted as a first approximation of the drag coefficient in open-channel flows.

Calculation of water surface profile for flow with emergent cylinders
Before examining the validity of Eq. (18) for calculating water depth profiles in gradually varied flows in conjunction with Eq. (6), we examine the accuracy of the equations for the uniform subcritical flows (L1-L5), as visualized in Fig. 11. Although the coefficient β Cp varies with respect to the vegetation density, as indicated in Fig. 7, the coefficient is set to β Cp = −0.4 in the upcoming analysis because the variation with respect to the density does not appear to be significant when the vegetation density is greater than 0.04. Furthermore, we want to focus on the effects of the non-equilibrium open-channel flows at the constant-density condition λ = 0.067. The water surface profiles are calculated using Eq. (6) from the downstream end water depth at x = − 0.1 m. The dragforce coefficient C D is calculated with Eq. (18) for Run U1, whereas the constant C D = C D0 is applied for Run U2. The RMSEs between the results of the calculation and experiment for L3-L5 are larger than those for L1-L2, indicating that the RMSEs depend on the velocity head. For the calculations, C D0 = 1.22 is adopted, which is obtained earlier in the study for Fig. 8 Fig. 11 Calculation results based on Eq. (6) and RMSEs of water depth for water surface profiles in uniform-flow conditions. Equation (18)  assumption introduced for Eq. (6). The velocity head causes the local water level to rise and is considered to be the representative length of the water surface measurement error, measured manually using a point gauge. The RMSEs are also affected by errors in the calculation model based on Eq. (6) from the momentum correlation coefficient in the advection term in Eq. (1) induced by the spatiotemporal variation in the velocity. The RMSEs of Run U1 are smaller than those of Run U2 except for Case L3, in which the effect of Eq. (18) for f s on C D is manifested in cases with relatively large F . The experimental data and discussions in this paper regarding water depth are considered to include an order of error similar to the velocity head, i.e., a few millimeters at the maximum. The water surface profiles for L1-L5 are well reproduced within this scale. Subsequently, the same calculation method with the same coefficients as those in Run U1 and U2 were applied to gradually varied flows, including accelerating and decelerating flows (Fig. 12). The calculation results from Eq. (6), in conjunction with Eq. (18) for C D and constant C D = C D0 , are indicated by dashed lines and small square symbols, respectively. Whereas the water depth profiles for the decelerating flows (b) were properly reproduced by the calculation results, those for the accelerating flows (a) were underestimated. Although the calculation method was designed to consider an error of a few millimeters even for uniform-flow conditions (Fig. 11), the underestimations by the calculations (Fig. 12a) were still significant relative to the error. On the other hand, compared to the differences between the measurement and calculation results of Eq. (6), the differences in the calculated results between Eq. (18) for C D and constant C D = C D0 were negligible compared to the deference from the measured results because relatively high Froude number flow restricted to near the boundary. This indicates that accelerating flows involve a mechanism that increases the drag coefficient or resistance of the obstacles. As shown in Table 2, the pressure measurement around a cylinder also revealed that the base drag coefficient C D0 for the accelerating flow (C D0 = 2.68) was larger than that for the uniform flow (C D0 = 1.79). The base drag coefficients are calculated using the pressure coefficient C p , which excludes the effect of the water surface from the coefficient, even though the coefficients β Cp for accelerating flow include the effect of the pressure gradient, as mentioned in a previous section. Figure 13 compares the distributions of the pressure coefficients for a solitary cylinder in a free stream [30], and for a single cylinder in a cylinder group in uniform and accelerating open-channel flows. The C p values of these experiments were calculated based on an assumption of C p = 1 at the stagnation point [44] to obtain an undisturbed water surface level, which is necessary for calculating the Δh distribution. Whereas the distributions of the positive pressure coefficients in front of the cylinder among the three cases are similar, the negative pressure coefficients in the separation zone differ. For the accelerating flow, the pressure in the separation zone decreases to values lower than those for the uniform flow, contributing to an increase in the drag coefficient for the accelerating flow. The mechanism whereby the pressure decreases in the separation zone is explained in Fig. 14, which compares the longitudinal water surface profiles for uniform and accelerating flows. It is assumed that, in the separation zone, the velocity is low [12,39,40], and the water depth is nearly constant; furthermore, there is exposure to the water depth just downstream from the separation zone. Therefore, the pressure in the separation zone is less for the accelerating flow than for the uniform flow because of the pressure gradient in the direction of flow.
The pressure reduction in the separation zone investigated in the experiment and visualized in Figs. 13 and 14 indicates that the additional drag force associated with the pressure gradient, f p , in the accelerating flow should be considered in evaluating the drag force f d of Eq. (8). The additional pressure due to the decrease in the separation area is expressed as where L is the representative length of the separation zone. Equation (25) is not attributed to the longitudinal buoyancy, which is calculated using dL = πd 2 /4, but indicates an increasing form drag associated with a decrease in pressure in the separation zone for the negative pressure gradient. Based on an assumption that the representative length is expressed by the cylinder dowel diameter d as L = kd (k, a coefficient), and a substitution of Eq. (25) into Eq. (8), the momentum equation for Eq. (1) can be re-written for steady flow conditions as follows: Equation (26) with a constant unit-width discharge (Eq. 2) yields the following depth profile equation for a flow with emergent vegetation: To calculate accelerating flows using Eq. (27), the coefficient k should be determined. In this study, the value k was obtained to minimize the discrepancies of water depth between the measurements and calculations based on Eq. (27). Figure 15 shows the RMSEs of the calculation results at all measurement points against k for all accelerating and decelerating flows. For all cases, the RMSE decreases from k = 0 and then increases after reaching a minimum, except for L10 in decelerating flow, which takes the minimum value at k = 0. The k value yielding the minimum RMSE varied by case, ranging from 0.0 to 2.1. Therefore, for each case, the k value that yielded the minimum RMSE is applied to the calculation using Eq. (27). According to Fig. 15, the calculations based on Eq. (27) with appropriate values k yielded low RMSEs in all cases, and the minimized RMSEs, i.e., 0.29-1.23 mm, were not significant relative to those in uniform flows (Fig. 11). The water surface profiles calculated using Eq. (27) with appropriate k values from Fig. 15 were then compared with those based on the measurement results, as shown in Fig. 12. The results obtained using this calculation method are in good agreement with the measurement results. Although introducing Eq. (25), for the calculation of f p , to the examination of decelerating flows had little effect on the validity of the analysis, the introduction of the additional drag force f p from Eq. (8) into the examination of accelerating flows was demonstrated to be a necessity. It must be recalled that the calculations with Eq. (6), which adjusts the drag coefficient for the uniform flow and accounts for water surface effects with Eq. (18), could not reproduce the water depth profiles for non-uniform flows. Because the practical use of the numerical calculation of flows with obstacles necessitates evaluation of the drag force for non-uniform flows, it is significant that several non-uniform flows can be reproduced using a unified calculation method. It is known that the relative separation length to cylinder diameter varies depending on the hydraulic conditions, including the Reynolds number and obstacle arrangement condition [12,39,40], whereas the representative length introduced to evaluate the additional drag force in Eq. (25) is indirectly related to the separation length. On the other hand, with regard to the proper evaluation of L with the appropriate value of k, that remains a question for future research efforts to resolve.

Conclusions
This study proposed additional terms for the drag-force expression with the drag coefficient to include the effects of the water surface and the pressure gradient in an average scale. The effects of adding these terms were investigated through two types of laboratory experiments on uniform and non-uniform open-channel flows with emergent cylinders, focusing on the entire longitudinal water depth profile and the distribution of local pressure around the cylinder. The main results of the study are as follows: (1) Based on the assumption of hydrostatic pressure distribution, it was determined analytically that the drag coefficient C D decreases in proportion to the square of the Froude number. The proportionality factor in the derived formula is composed of the drag coefficient of the base component of the drag force, C D0 , and the second-order moment of the pressure coefficient, β Cp . The derived formula can reproduce the characteristics  of the drag force on emergent cylinders in a uniform flow, as obtained from our experiments and from others previously conducted. (2) For low Froude number flow, the water surface effects on the drag force for the depth profile calculation were negligible in this study. Although the drag coefficient formula for uniform flows can be applied to the calculation of water depth profiles in decelerating flows with emergent cylinders, the same formula led to an underestimation of resistance in accelerating flows, where the pressure in the separation zone decreased because of the negative gradient of the water depth in the longitudinal direction. Therefore, a modified method with additional drag-force terms accounting for the pressure drop in the separation zone was proposed. The accuracy of the acceleration flow reproduced therewith was improved considerably through the incorporation of this additional drag force into the calculation.
As reported in several past publications on the drag force in open-channel flows, the drag coefficient depends on hydraulic conditions and on several factors pertaining to the obstacles, such as their arrangements and shapes. Thus, estimating an appropriate value for the drag coefficient of the base component of the drag force is challenging. In fact, parameterization of the channel resistance coefficients is inevitably necessary in practical applications such as flood flow simulations of rivers [19,20]. The most practical use of numerical calculations for flows with obstacles necessitates an evaluation of drag forces in non-uniform flows. The main contribution of this research is the development of an appropriate modification of the drag force from the base component for several non-uniform-flow conditions using a simple unified expression, which accounts for the open-channel effects due to water surface variation and pressure gradients in accelerating flows on the drag-force coefficient. It is expected to be particularly useful for resistance evaluation of obstacles in high Froude number flows with strong non-equilibrium and unsteadiness, such as dam-break flows and tsunamis, where both effects of the water surface and pressure gradient on the drag force are considered to be significant. In order to improve the validity of the analysis of the nonequilibrium flow field, it is necessary to clarify the effect of the pressure gradient on the drag force, i.e., the representative length L.