Transverse flow distributions in ice‐covered channels

Traditional discharge measurements in ice‐covered river channels are time‐consuming and costly, thus limiting the elucidation of flow patterns during winter. This study proposes an analytical model for quantifying the transverse distributions of the unit‐width discharge and the depth‐averaged streamwise velocity in ice‐covered channels by extending an existing stream‐tube method to include the effects of cover and bed roughness and cross‐sectional geometry. A new set of equations are obtained for this extended method. An empirical coefficient α is included to lump the effects of the shape and friction factors of the ice cover and the channel bed. The model is compared with flume experiments and field data for fully and partially ice‐covered flows to examine the variation of the coefficient α.


| INTRODUCTION
The transverse flow distribution under ice-covered conditions could affect the transport and mixing of pollutants, suspended materials, sediment transport, and morphological changes in river channels (Sayre & Song, 1979;Shen, 1978;Yotsukura & Sayre, 1976). Additionally, it provides information on the river thermal-ice processes, such as the evolution of undercover transport of frazil granules that contribute to the development of hanging dams (Shen & Wang, 1995) and the development of border ice (Huang et al., 2012). The existence of ice cover provides thermal insulation and additional surface resistance for river flow, and thus significantly changing the boundary shear stress distribution, velocity profiles, and turbulence characteristics (Easa, 2017;Koyuncu & Le, 2022;Tatinclaux & Gogus, 1983). Teal et al. (1994) systematically compared 10 pointvelocity measurement techniques with field velocity profiles under the ice cover. They suggested that velocities at relative flow depths of 0.2 and 0.8 can provide good estimation of the depth-averaged velocity and developed the two power law for vertical velocity distribution by treating the ice-covered flow as two open channel flow. Demers et al. (2013) compared velocity profiles from open channel and icecovered conditions collected in Mitis River, Canada. They found that the ice cover and the river bed generate and sustains two sets of macroturbulent coherent structures which compete with each other in the lower and upper layer and affecting the vertical velocity profiles.
Although ice covers significantly impact river flows, studies on transverse flow distribution in ice-covered rivers are limited due to the challenging field conditions and timeconsuming operation to obtain good velocity data. Shen and Ackermann (1980) proposed a stream-tube method to estimate the transverse unit-width discharge distribution in rivers with different ice cover conditions and channel characteristics. Yang (2015) developed an analytical model to estimate the depth-averaged streamwise velocity based on the Shiono and Knight Method (SKM) for open channels (Shiono & Knight, 1991). The SKM has been further applied to partially and fully ice-covered flows in rectangular and compound experimental flumes (F. Wang et al., 2020;Zhong et al., 2019). Their predicted results are reasonably consistent with laboratory flume experiments. However, ice and bed friction factors are required in addition to the cross-sectional geometry to apply SKM. In this study, the stream-tube method of Shen and Ackermann (1980) is examined with additional laboratory and field data for possible improvements.

| THEORETICAL ANALYSIS
In this section, the flow distribution in ice-covered will be analyzed considering cover and bed roughness, which was ignored by Shen and Ackermann (1980). Uzuner (1975) reviewed methods for calculating the velocity and the composite bed and ice cover flow resistance in channels and discussed the validity of the composite Manning's coefficient and the logarithmic velocity profiles (P. Larsen, 1973;P. A. Larsen, 1969). For analyzing flow in ice-covered channels, the ice-covered channel cross-section is divided into the ice and the bed zones, as shown in Figure 1. The theoretical division line is the line of zero shear stress at the interface between the ice and bed influencing zones. The zero-shear stress line is considered identical to the maximum velocity line as their differences are negligible for engineering applications in natural rivers (Krishnappan, 1984;F. Wang et al., 2020;Zhong et al., 2019). For ice-covered streams, the total discharge is then divided into two parts, Q i and Q b , representing the cover and bed affected portions of the discharge, as:

| Governing equations
where Q is the discharge, f is the Darcy-Weisbach friction factor, A is the flow area, = R A P , is the hydraulic radius, P is the wetted perimeter, S is the energy gradient, g is the gravitational acceleration, and subscripts i and b denote the hydraulic parameters in the ice and the bed zones, respectively.
Combining the flow in the ice and bed zones, the total discharge under the ice cover becomes: Similarly, the discharge through the partial cross-section , where the subscript p denotes variables in the partial cross sectional area. The nondimensionalized partial discharge becomes where β is the ratio of the ice cover and river bed friction factors. The derivation of Equation (3) is detailed in the Appendix, which shows that Q Q p is affected by the ratios of flow area, hydraulic radius, wetted perimeter, and the ice and bed friction coefficients. Equation (3) can be rewritten by introducing an empirical coefficient α as: in which the coefficient α is: Physically, the coefficient α is affected by the crosssection geometry and the ratio of ice and bed roughness. Shen and Ackermann (1980) obtained a constant α = 2/3 by assuming the mean velocity in A i and A b are approximately the same, and n i ≃ n b (n is noted as Manning's roughness coefficient) or P i ≃ P b . The calibration and sensitivity analysis of α will be discussed in Sections 3 and 4 based on laboratory and field data. Equation (5) might be used to estimate the value of α with known ice and bed friction factors. However, these factors, especially the ice cover friction factor, can vary during the winter (Ghareh Aghaji Zare et al., 2016;Shen & Yapa, 1986). It also shows that when the ice and bed friction factors are the same, the coefficients β = 1, and α= 1/2, independent of the cross-sectional shape. Additionally, since the interfacial shear is neglected in the above derivation, it is necessary to use the average of Q p calculated from both left and right banks in the flow distribution analysis.
With the dimensionless partial discharge distribution in Equation (4), the transverse distributions of the unit-width discharge along the cross-section can be calculated through numerical differentiation. The unit-width discharge at a given transverse location, y, along the cross-section is: and the dimensionless unit-width discharge is: where, ̅ q y is the dimensionless unit-width discharge, and B is the river width. Substituting Equation (7) into (6) yields: Using the measured unit-width discharge, q , y and calculated ̅ q y , Equation (7) can be used to estimate the total flow discharge. The depth-averaged streamwise velocity is then rewritten as: where h y ( ) is the flow depth at the location y. Equations (7), (8), and (9) are the main equations for estimating discharge and streamwise velocity distributions in both open and ice-covered flows, similar to the procedure given in Shen and Ackermann (1980), but consider the effect of α coefficient. The cross-sectional geometry, that is, the flow depth distribution can be obtained either by depth measurements through drilled holes on the cover along the cross-section or using the double-frequency ground-penetrating radar (Fu et al., 2018). Equation (7) can be used to calculate the total discharge with one measured unit-width discharge. Equations (8) and (9) provide the dimensionless unit-width discharge and the depth-averaged streamwise velocity distributions with known flow depth distribution.
In the following analysis, the vertical velocity profiles along the thalweg is used to calculate q y0 and three additional velocity profiles at the left berm, the main channel, and the right berm are used to calculate three α values, respectively. The average of these three values provides the calibrated α at a specific cross-section.

| VALIDATIONS WITH FLUME EXPERIMENTS
Although flume experiments may not reflect the complexity of field conditions, they provide accurate data sets for validating analytical formulations. This section uses a set PAN ET AL.

| 199
of flume experiments conducted at China Institute of Water Resources and Hydropower Research (IWHR) and F. Wang et al. (2020) to validate the analytical formulation.

| Experimental setup
A series of experiments were conducted in a 17.3 m long tilting Plexiglas flume of IWHR Hydraulics Laboratory, with a compound cross-section, as shown in Figures 2 and 3. The ice cover was replicated with 16-m-long and 1.18-m-wide smooth foamed plastic. All experiments were conducted at room temperatures of about 20°C, considering that the thermal effect is negligible. The flow was recirculated using a centrifugal pump and a mechanical valve. A flow straightener was used to straighten the flow at the inlet. A DN100 KRONE electromagnetic flowmeter was used to provide steady flow discharge with the upper limit of 0.11 m 3 /s. The flowmeter accuracy is ±0.5%. The flow depth was controlled by a tailgate shown in Figure 3a. Seven scales were attached to the sidewall of the flume to record water levels at a longitudinal interval of 0.25 m. Figure 3b shows the cross-sectional shape of the channel with a rectangular main channel and two symmetrical berms. The bed elevation difference is 0.075 m. There are 23 velocity gauging locations spaced at 0.01, 0.05, or 0.10 m intervals to measure vertical velocity profiles across the channel for transverse flow distributions in each test run. Smaller intervals were used near the wall or the transition between the main and the berms to capture velocity variations accurately.
Vertical velocity profiles were measured at the crosssection 12 m from the cover leading edge, where the uniform flow was fully developed. These point velocities were measured with a SonTek 3D Acoustic Doppler Velocimeter (ADV) at a vertical interval of 0.01 m. The standard measurement accuracy of the ADV is ±0.5%, while the vertical distances measurement accuracy is ±0.001m. The ADV has a 0.01 m blind zone close to the river bed and a 0.02 m blind zone near the bottom of the ice cover. Experiments were conducted under uniform flow conditions with discharges varying from 0.0208 to 0.1040 m 3 /s, as listed in Table 1. Figure 4 shows two transverse discharge distributions from the measured velocity data. Since the flume wall and the simulated ice cover are both smooth with identical friction factors, the parameter β = 1. The parameter α is calibrated to be 0.5, which is consistent with Equation (5). The model reproduces unit-width discharge distributions for the compound channel with rapid depth changes. The average relative error is about 6.8%, mainly from the interfacial shear stress effects at the rapidly changing flow depth locations.
F. Wang et al. (2020) measured the depth-averaged streamwise velocity for a flume with a symmetrical trapezoidal cross-section with a sidewall slope of 45°for fully covered and symmetrical border ice conditions. Data obtained by F. Wang et al. (2020) was used to validate the analytical model with calibrated α = 0.25 for partially and fully covered conditions. Figure 5 shows the calculated results from the proposed model compared with the measured data and results from the SKM model proposed by F. Wang et al. (2020). The comparison shows half of the symmetrical trapezoidal cross-section starting from the centerline of the channel. The calculated depth-averaged velocity is consistent with the measured data and calculated from the SKM model of F. Wang et al. (2020) for both fully and partially covered conditions. The average relative error of the proposed model was approximately 3.4%. There are significant deviations between calculations made by F. Wang et al. (2020) and this study near the right wall, where sidewall effects become important. Unfortunately, there were no detailed measurements near the wall for comparison.

| APPLICATION TO FIELD CONDITIONS
The present analytical model is applied to field data from the Inner Mongolia reach of the Yellow River with full and partial ice cover conditions. To simplify the model's application, only one value is calibrated for α at a given cross section.
4.1 | Comparison with field data C. Wang et al. (2017) observed the unit-width discharge profiles in the winter of 2014-2015 at four cross-sections in the Inner Mongolia Reach of the Yellow River, along with cross-sectional geometries and ice-cover distributions at different times, as shown in Figures 6-9. The average discharge varied from 300 to 700 m 3 /s during the observation period, while the ice thickness changed from 0.3 to 0.6 m. The top surface of the ice cover was set as the reference level = z 0. The ice cover thickness and bed geometry were referred to as the z axis. Figures 6 to 9 compared the observed transverse flow distribution with the present method and that of Shen and Ackermann (1980). Figure 6 compares the calculated and measured unit-width discharge distributions for ice-covered flow at the Wenbuhao on February 15 and March 13, 2015. In this case, α is calibrated as 2/3 while the aspect ratio is approximately 20. Calculated unit-width discharges agree well with field data. The maximum unit-width discharge decreases by more than 20%, from 660 to 485 m 3 /s, from February 15 to March 13, as shown in Figure 6a,b. The bed elevation at the thalweg increased by about 2 m in a month, showing that sediment transport and bed elevation changes are significant even in fully ice-covered periods. The model works well under such conditions. Figure 7 shows the measured and calculated unit-width discharge, and the observed bed and ice thickness distributions at the gauge of Dabusutai on January 23 and February 9. With calibrated α = 0.5, the calculated unit-width discharge agrees well with the observation for a wide variation of flow depth from 0 to 9.5 m. Figure 8 shows the calculated results, observed unit-width discharge, bed, and ice cover thicknesses at Baotou on  January 30 and February 13 with α = 0.5. The present model performed well in estimating unit-width discharge distributions under complex cross-sectional shapes in fully ice-covered flows. The empirical coefficient α remains constant, 0.5, at Dabusutai and Baotou with changing ice cover thickness and bed geometry at different times. The model was further applied to partially ice-covered flows at Sanhuhekou with α = 0.5, as shown in Figure 9. The model again works well compared with measurements in nonsymmetrical border ice conditions. The ice cover and the bed elevation varied significantly from January 9 to January 19 during the freeze-up period at Sanhuhekou. The dynamic ice and alluvial bed verify the capability of the analytical model. There are some deviations in the calculated discharge at Sanhuhekou on January 9 at the transition zone from border ice to open flows. These deviations may decrease if α is calibrated separately for the ice-covered and open flow zones. Overall, the present model performs well in the Yellow River with complex and varying bed and ice cover profiles. The comparison of unitdischarge distributions with the Shen and Ackermann (1980) method showed significant improvements especially for cross sections with sharp changes in depth and partially ice-covered condition. The comparisons showed the present method could improve the maximum unit-width discharge by about 20%, as in the case of Figure 7. This is important when using the method to calculate the total discharge based on the maximum unit-width discharge.

| Sensitivity analysis
The unit-width discharge and the depth-averaged streamwise velocity distribution are functions of ice and bed friction factors, ice and bed influencing area, and crosssectional geometry in the winter ice season. Their transverse distributions are calculated through the improved stream tube method with empirical coefficient, α. This coefficient lumps all uncertainties from channel geometrical irregularity and nonuniform boundary roughness of ice and channel bed. The α coefficient was found to be essential for the present model. The calibrated value of α = 2/3 at Wenbuhao, which has an aspect ratio of 20, is the same as that of Shen and Ackermann (1980). However, α is calibrated to be 0.5 for the other three cross sections at Dabusutai, Baotou, and Sanhuhekou under both fully and partially ice-covered conditions. Their aspect ratios exceed 50, which is 1.5 times larger than that at Wenbuhao. The aspect ratio is approximately 5 in the two flume studies, and α varies from 0.25 to 0.5.

| 203
Further sensitivity analysis of the model on α was conducted in the IWHR flume experiment case 4 and the Yellow River case at Sanhuhekou. Figure 10 shows variations of estimated unit-width discharges with α values varying between 0.25 (50% decrement) and 0.67 (34% increment) from the 0.50 calibrated value. The unit-width discharge distribution becomes more uniform with decreasing discharge at the main channel and increasing discharge at berms when α decreases from 0.67 to 0.25. Figure 11 shows the sensitivities of the model to α in the Yellow River case at Sanhuhekou on January 19, 2015. The unit-width discharge at the thalweg decreased by 14% with a 50% decrease in α, while that of the shallow area decreased by 10% with a 34% increase in α.
As indicated by Equation (5), α varies with different ice and bed friction factors and cross-sectional shapes. Figure 12 shows the calculated α from Equation (5) at different transverse locations at Sanhuhekou on January 19 for β of 0.5 (50% decrement on the reference value), 1.0 (reference F I G U R E 9 Comparison of transverse unit-width discharge distributions for partially ice-covered flow at Sanhuhekou: (a) January 9 and (b) January 19, 2015.
F I G U R E 10 Sensitivities of estimated unit-width discharge to coefficient α in experimental case 4.
F I G U R E 12 Variations of α calculated via Equation (5) with different ice and bed friction factor ratios β at Sanhuhekou on January 19, 2015.
F I G U R E 11 Sensitivities of estimated unit-width discharge to coefficient α in the Yellow River at Sanhuhekou on January 19, 2015.
value for the same ice and bed friction factors), and 1.5 (50% increment on the reference value). Theoretically, the value of α is a constant of 0.5 when the ice friction factor is equal to the bed friction factor. Under this condition, α is independent of the geometrical irregularity as shown in Figure 12. Otherwise, α varied along the transverse axis y at Sanhuhekou for different flow depths and ice thicknesses, owing to shape effects. Figure 12 also shows that α increases with any increase in β in the border ice regions.
It is suggested that α should be calibrated with the measured unit-width discharge at three or more points for a cross-section as α remains constant at different times in Figures 6-9. The preliminary study shows that α varies from 0.25 to 0.67, with the typical value of 0.5 for general applications. However, the range of α might change in natural rivers considering the change of ice cover roughness over the season (Shen & Yapa, 1986). Ghareh Aghaji Zare et al. (2016) showed that the ratio of ice to bed roughness varies from 0 to 4 during the 4-month winter period in Lower Nelson River, Canada. For channels with known bed roughness, the cover roughness could be obtained concurrently with the discharge measurement by observing the stages at two cross sections at a short distance up and downstream of the discharge measurement cross section.

| CONCLUSIONS
In this study, an analytical formulation for the stream tube method on transverse flow distribution with ice effects is developed to examine the effect of ice cover roughness and cross-sectional geometry and applied to data from laboratory  flumes and natural rivers with aspect ratios varying from 5 to 50. The formulation is based on the two-layer flow partition and the Darcy-Weisbach friction factors of the ice cover and river bed. The formulation include a coefficient α related to the cross-sectional geometry and the ratio of ice and bed roughness. The model reduced to the Shen and Ackermann's (1980) method when α is 2/3. The model uses only one vertical velocity profile for calculating the total discharge, unit-width discharge, and depth-averaged streamwise velocity in open, partially ice-covered, and fully ice-covered flows. However, coefficient α remains empirical in this study owing to the limited availability of databases for ice-affected flow. Further studies are required to elucidate α in terms of flow variations, river geometries, and boundary roughness.  (1980)'s stream tube method, the interfacial shear stress is assumed to be negligible, and the secondary flow has no net influence on the depth-averaged streamwise velocity. For steady uniform flow, the gravity component along the channel slope is balanced by the shear forces from the ice cover and the river bed.
where τ is the mean boundary shear stress, and ρ is the fluid density. Equation (A1) implicitly assumes that the river slope is mild. As the partial flow area is divided into the ice cover and river bed dominating zones, and there is no shear stress along the theoretical division line, the force balance equations for these two zones are as follows: The total energy slope is assumed to be the same as the energy slope in the ice and bed layers: The ratio of the ice cover and river bed friction factors is defined as: Furthermore, the mean velocity is assumed to be the same in ice and bed zones (Uzuner, 1975). Therefore, the ratio of the ice-covered shear stress and the bed shear stress is as follows: where U is the cross-sectionally averaged velocity. Combining Equation (A2) with (A6) yields the following equation: Equation (A7) gives: Dividing the discharge through the partial cross-section area with total discharge yields: