Numerical simulation of mechanical breakup of river ice-cover*

WANG Jun (王军), HE Liang (何亮), CHEN Pang-pang (陈胖胖)

School of Civil Engineering, Hefei University of Technology, Hefei 230009, China, E-mail: junwanghfut@126.com

SUI Jueyi

Environmental Engineering Program, University of Northern British Columbia, Prince George, Bristish Columbia, Canada

Numerical simulation of mechanical breakup of river ice-cover*

WANG Jun (王军), HE Liang (何亮), CHEN Pang-pang (陈胖胖)

School of Civil Engineering, Hefei University of Technology, Hefei 230009, China, E-mail: junwanghfut@126.com

SUI Jueyi

Environmental Engineering Program, University of Northern British Columbia, Prince George, Bristish Columbia, Canada

(Received January 26, 2013, Revised March 31, 2013)

Ice jams and ice dams in rivers will cause significant rises of water levels. Under extreme conditions, the ice flooding during winter or early spring may occur. In this paper, by considering the fluid-solid coupling effect caused by the water and the ice cover, the mechanisms of the mechanical breakup of the river ice cover are studied. A formula is obtained for determining whether or not the mechanical breakup process would happen under the hydraulic pressure of the flow. Combined with the hydraulic model under the ice covered flow, a numerical model is built and the interaction between the discharge, the hydraulic pressure under the ice cover and the date for the mechanical breakup of the river ice cover is simulated. The simulated results of the dates for the mechanical breakup of the river ice cover agree very well with the field observations of the breakups of the river ice cover in the Hequ Reach of the Yellow River. Therefore, the numerical model might serve as a good preliminary step in studying the breakup of the river ice-cover, evidencing many important parameters that affect the ice-cover process.

criterion for breakup, fluid-solid coupling effect, Hequ Reach of the Yellow River, ice cover, mechanical breakup, numerical model

Introduction

In winter, the river ice represents an important hydrologic element in temperate and polar environments. An ice cover alters the hydraulics of an open channel by imposing an extra boundary on the flow, altering the water level and the flow velocity as compared to those under the ice free conditions. A major consequence of the ice-cover formation in many rivers is the associated jamming, particularly, during the periods of the fall freeze-up or the spring break-up. In fact, the ice jamming is one of the most conspicuous and momentous ice-related phenomena. Ice jams can result in significant rises of water levels. Moreover, the backwater effect of an ice jam on the water level increases with the jam thickness. The increases in the jam thickness, in turn, increase the channel slope and the flow velocity, and may breakup the river ice cover. When the leading edge of an ice cover or a jam progresses to the section with a high flow velocity, the incoming ice floes will (1) either submerge and become entrained on the underside of the cover or (2) result in the mechanical breakup of the river ice cover in the downstream river reach.

Ice dams may be created during the thawing periods of the spring months. If certain parts of a waterway remain frozen while the rest is thawed, the lingering ice can create a barrier. In spring, the ice cover in some rivers that flow north starts to breakup from upstream to downstream due to the difference of the temperatures between the upstream and the downstream (such as the Lower Reach of the Yellow River in China). The breakup of the river ice cover in the upstream river reach results in an increase of the flow discharge and a rise of the water level in the rivers. Since the ice cover in the downstream reach of the rivers (north-flowing rivers) is still continuous, thehigh kinetic energy of the flowing water and the ice floes due to the breakup of the upstream river reach may result in the mechanical breakup of the river ice cover in the downstream river reach. Generally, the mechanical breakup of the river ice cover will produce a vast amount of ice floes. Depending on the characteristics of the downstream river reach, the ice floes will be blocked to form ice jams and hanging dams in the downstream river reach.

The most common location for the formation of ice jams and hanging dams is places with reduced or non-uniform velocities. Obvious locations include where the river enters a lake or reservoir, along a shallow water river reach, along a braided river reach, in front of the leading edge of the ice cover, at river bends or when the flow encounters a deep pool. The formation of jams can be particularly damaging when the previously mentioned features are combined with an upstream stretch of rapids where large quantities of frazil are produced. When an ice dam breaks, it can flood the areas downstream causing damages in the surrounding natural and human environments.

Due to their inherent nature, the ice problems are invariably complicated. The following factors must be considered in the studies of river ice problems:

(1) River morphology, including the channel slope, the channel geometry, the channel networks, the shoals, the riverbed roughness, the vegetation in the channel and the floodplain, and others.

(2) Hydraulic conditions, including the flow velocities, the water depth, the water surface profile, and the downstream conditions.

(3) Meteorological factors, including the air temperature, the water temperature, the wind speed, and the snowfall.

(4) Many other factors, including the human interventions, the ground heat input, and others.

In the past three decades, the variations of water levels under ice-covered conditions[1-5]and the evolutions of ice jams[1,2,6]were amply studied. On one hand, due to a large aggregate ice thickness and a high hydraulic resistance relative to the sheet ice, the ice jams tend to disturb the riverbed[7,8], on the other hand, the ice jams can cause mechanical breakup of the river ice cover.

So far, the research work was focused on the mechanical breakup of the river ice cover[9-11], the impact of the thermodynamics on the ice cover was not duly considered, and no adequate numerical models were developed. To forecast the river breakup, either the methodology of the neural network[12,13]or the empirical models[14]were employed. The objective of this paper is to derive a formula for determining whether or not the mechanical breakup process would happen under the current hydraulic pressure. Combined with the hydraulic model for the ice covered flow, a composite model is developed to simulate the interaction between the discharge, the hydraulic pressure under ice cover and the date of the mechanical breakup of the river ice cover.

1. Numerical model

1.1 Models for criterion for onset of mechanical

breakup

The breakup of an ice cover often results in the formation of ice jams. However, there are not many reliable methods for predicting the ice cover breakup. One method uses the hinge resistance at the riverbanks as a criterion for the onset of the ice cover breakup in relatively straight rivers[15]. The shear stress on the ice cover can be estimated with reasonable accuracy[15], and then the counteracting hinge resistance might be determined at the onset of the cover breakup. However, the calculated value of the hinge resistance does not vary much year in and year out, and thus the proposed criterion appears to be only applicable for the breakup of the ice cover in other river reaches[15].

In winter or spring, the changes of hydraulic conditions in an ice-covered river lead to the changes in the pressure force acting on the ice cover. As a result of the increase in the hydraulic pressure acting on the ice cover due to the rising water level, a stress concentration will be created along the edge of the ice cover near the river banks, which will lead to the formation of cracks (longitudinal and transversal cracks) in the ice cover and the river ice cover may break. The formation of longitudinal cracks in an ice cover is a sign that the river ice cover may break. Thus, in this study, the formation of a longitudinal crack in an ice cover will be used as the criterion for the mechanical breakup of an ice cover.

In general, the flow in a straight channel is normally a symmetric flow, if we use the longitudinal centre line of the channel as the axis, it will be an axis-mirror symmetric flow. Conceptually, under the ice covered conditions, the force acting on the ice cover should also obey the “axis- mirror symmetry”.

Fig.1 Schematic diagram of load on the ice cover and its deflection

Figure 1 shows the schematic diagram of the load acting on an ice cover and the deflection of the icecover from the centre of the channel to the riverbank. The load includes the following important forces: the longitudinal drag force from the flowing water (in the stream-wise direction), the uplift force from the water (in the vertical direction) and the weight of the ice cover. As pointed out in Ref.[15], since the longitudinal crack in an ice cover is mainly caused by the vertical uplift force, the effect of the longitudinal drag force on the formation of the longitudinal cracks will not be considered in this study. Additionally, the forces in the flow direction are perpendicular to the longitudinal crack. Thus, the hydraulic factors such as the flow velocity and the hydraulic slope need not be considered so far as knowing the stress force on the bottom surface of the ice cover in the hydrodynamic model.

Before the formation of longitudinal cracks in an ice cover, the ice cover is fixed to the riverbanks[9,10]. With the increase of the discharge in the channel, the water level rises. As a result, the distributions of the uplift force(p)in both the longitudinal and transversal directions will be changed. Compared with the change in the transversal direction, the change in the longitudinal direction can be ignored. Similar to the assumption in Ref.[11], in this study, the uplift force (p)in the transversal direction is assumed to be distributed uniformly along the studied ice cover. As indicated in Fig.1, the AB represents the ice cover, the x-axis is along the transversal direction and pointing to the centre of the channel, Point A represents the point near the river bank, Point B represents the centre of the ice cover in the transversal direction (note: near the centre of the channel, there is no further deflection of the ice cover since all forces at this point are in equilibrium),qArepresents the load at Point A, and l represents the horizontal length of the studied ice cover AB.

In this study, the deformation of the ice cover is considered. The forces acting on the bottom of the ice cover include the uplift force caused by the water (such as the rises of the water level) and the buoyancy force, as well as the weight of the ice cover. The sum of these forces per unit area (pT)can be described by the following equation,

where γand γiare the specific weight of the water and the ice, respectively,hiis the thickness of the ice cover before the breakup,y is the deflection of the ice cover, and siis the relative mass density of the ice, si=γi/γ.

Equation (1) can be simplified as Along the transversal direction, the load(q)acting on the ice cover is due to the total uplift force (pT). Here, the load is assumed to be uniformly distributed along the transversal direction as shown in Fig.1, which means that the load from the uplift force(p)on the bottom of the ice cover is uniformly distributed. The buoyancy term,γ(sihi-y), decreases with the increase in the deflection of the ice cover in the vertical direction. The load caused by the buoyancy term is a gradually varying load. Since the thickness of the ice cover in a straight river reach does not change much along the transversal direction[1,2], the load from the weight of an ice cover can be considered uniformly distributed in the transversal direction.

Fig.2 Analysis of force on the studied ice-cover element

As shown in Fig.1, the total load(q)from all individual forces mentioned above decreases alongxdirection (from the river bank to the centre of the channel). As shown in Fig.1, the distribution of the total load(q)in thex-direction (the transversal direction) can be assumed to be expressed by a linear function as

As shown in Fig.2, along thex-direction, the studied elementdx has a width of b(along the longitudinal direction).FQrepresents the shear force on the right side of the studied element. Since the studied element is in equilibrium, the external forceqand the internal force FQmust be in equilibrium. The following equation can be used to describe the shear force,

Since d FQ=qd x , from Eqs.(2)-(4), we have,In Eq.(5), if Point A is infinitely near the solid river bank, one can say that xA=0 and yA=0. Thus, Eq.(5) can be written as,

For a clamped-end beam[16], the differential equations for the deflection of the beam are as follows:

FQ(x)and M( x)represent the distribution of the shear force and the bending moment in thex-direction, where

The boundary conditions for a clamped-end beam can be expressed as[16]:yx=0=0 and y′x=0=0. Thus, the equation for the deflection curve of ice cover can be written as

In order to determine the deflection at B(yB),xis replaced byl in Eq.(7), then the deflection of the ice cover can be calculated as

whereI is the moment of inertia of the beam cross section about the neutral axis, and I= bh3/12. Ati Point B,PT=0. Then, the deflection of the ice cover at Point B is the maximum deflection of the ice cover in the transversal direction(ym).

From Eq.(2),p=γym, and from Eqs.(6) and (8), one has

From the differential equations for the bending moment and deflection curves, the maximum bending moment can be expressed as

In addition, the maximum bending strength of the ice cover should satisfy the following relation,

where W= bh2/6is the bending coefficient of azi cross section of the ice beam. From Eqs.(6), (9)-(11), we have,

where E0and σ0is the nominal value for bothE and σ, respectively,SσandSErepresent the source terms forEandσ, respectively. From Eqs.(12)-(14), we have,

Equation (15) can be used to calculate the uplift force during the breakup period. Conceptually, there should be one criterion for the river ice cover to break. During the breakup processes, since the discharge normally increases, the uplift force will increase. If the uplift force is above this criterion value, cracks along the river banks will appear. The formation of these cracks destroys the integrity of the ice cover with the river bank. The mechanical breakup of the river ice cover may result.

1.2 Hydraulic model for river ice cover breakup

In practice, it is difficult to determine the force on the bottom surface of the ice cover, since the force depends on various hydraulic variables. Also, the flow in natural channels is normally turbulent unsteady flow. The Navier-Stockes Equations can be used to describe the flow in natural channels. The governing equations are normally expressed as thek-εtwo equations module,

Equation (16) is the continuity equation, Eq.(17) is the momentum equation, Eq.(18) is thek equation, and Eq.(19) is theεequation, where uiis the velocity component inxidirection (i =1,2,3),ρis the volume-fraction based the average mass density,P is the corrected pressure,µtis the volume-fraction based average coefficient of the molecular viscosity, µtis the turbulent viscosity µt=ρCµk2/ε,Cµis an empirical constant of 0.09,σkand σεare the turbulent Prandtl numbers forkandε,σk=1.0 and σε=1.33, respectively,C1εand C2εare constants of theεequation,C1ε=1.44 and C2ε=1.92, respectively,Gkis the turbulent kinetic energy caused by the average velocity gradient,

In our previous research[4], the k-εtwo-equation turbulence model was used to simulate the ice accumulation under the ice cover. All simulated thickness values of the ice accumulation agree reasonably well with the measured values in laboratory. In Ref.[4], this hydrodynamic model was used to determine the lifting force on the bottom surface of the ice cover caused by current. The main advantage of the k-εmodel is that the k-εmodel can be used more precisely describe the impacts of the turbulent flow on the changes under hydraulic conditions. Also, to apply this model during the breakup of the ice cover, it is only necessary to assume that the borders of the ice cover are attached to the river banks. The hydrograph at the upstream entrance cross section of the river reach must be known. Also, both the thickness of the ice cover and the water level at the downstream cross section should be given. Then, the developed model will be solved to determine the stress on the bottom surface of the ice cover. Thereafter, one can determine (together with Eq.(15)) if the river breakup is a mechanical breakup.

Fig.3 Hequ Reach of the Yellow River

2. Use of the developed model for mechanical breakup of river ice cover

By conducting the regression analyses, the coefficient C in Eq.(15) is determined as C=2.0. Based on experimental results for ice in laboratory, Sui et al. (1994) obtained the following results by using a regressional method,σ0=450.61,Sσ=-3 .02Ti-86.73,E0=0.12and SE=-0.25. During the river breakup period, based on the field measurements at the Hequ Reach of the Yellow River as shown in Fig.3, the temperature of the ice cover is assumed to be –5oC in this study. A 3-D river model is developed based on a river bend section of the Hequ Reach of the Yellow River. The 3-D river model is divided by structured grids[18]. The pressure-velocity coupling is solved by using the SIMPLEC method. The second order upwind difference scheme is used to determine the average interpolated value for the kinetic energy, the turbulent energy and the energy dissipation rate. The interpolated value of the pressure fields is determined by using the BFW scheme.

Table 1 Forecasted dates of river breakup compared to the observed dates of river breakup

As for the boundary conditions used in this simulation, the following assumptions are made: (1) the standard wall function is used for describing the river bed and the bottom of the ice cover since both the riverbed and the ice cover are solid, (2) at the entrance of the cross section, the velocity of the flow is used as the entrance boundary condition, and (3) the pressure determined by the water level is used as the boundary condition at the exit cross section.

During the period before the river breakup at the Hequ reach of the Yellow River, with the increase of the temperature during early spring, the discharge in the Yellow River does not normally change much, and the fluctuation of the water level is very small. During the river breakup period, however, the discharge increases dramatically which results in a significant rise of the water level[2], which would lead to increases of the hydraulic pressure under the ice cover. From Eq.(2), if the deflection of the ice cover is zero, the

uplift force (p)is the upward resultant force caused by the hydraulic pressure force under the ice cover and the weight of the ice cover. This uplift force depends on the linear distribution of the hydraulic pressure on the bottom of ice cover. Therefore, one can use the increment of the pressure under the ice cover to represent the increase of the uplift force.

Using the data measured at the Shiyaobu cross section of the Hequ Reach of the Yellow River, the changes of the hydraulic pressure on the ice cover for different discharges are determined by combining the model for the mechanical breakup of the ice cover with the model for the river hydraulics under ice covered conditions. The calculated changes of the hydraulic pressure on the ice cover are used to determine the daily uplift force (p). By comparing the calculated uplift force(p)with the criterion for the mechanical breakup of the river ice cover, one can determine the date of the mechanical breakup of a river ice cover.

Table 1 lists both the calculated date of the mechanical breakup of an ice cover and the observed date of the mechanical breakup of an ice cover at the Hequ Reach of the Yellow River. It can be seen from this table that the calculated date for the river breakup in 1993 and 1994 is one day earlier than the observed dates. The following facts may be responsible for this difference: during the river breakup periods in 1993 and 1994, similar to other years, the discharge is gradually increased due to the increases of the temperature. Eight days after the temperature is above zero degree, the discharge is increased from approximately 400 m3/s to about 700 m3/s. However, the thickness of the ice cover in 1993 and 1994 is much larger than that in other years. As a result, when compared to the calculated dates for the river breakup by using the developed model, the observed dates for the mechanical breakup of this thicker ice cover in 1993 and 1994 is delayed by one-day. In overall, as long as the hydraulic conditions for the mechanical breakup of a river ice cover are reached, one can easily use the criterion model for the river breakup to determine the date of the river breakup.

3. Conclusion

In the present study, by considering the fluidsolid coupling effect caused by the water and the ice cover, the mechanisms of the mechanical breakup of the river ice cover are explored. The critical uplift force for the mechanical breakup of the river ice cover is identified. A model is developed as the criterion for determining whether or not the breakup process of the river ice cover is the mechanical breakup. By using laboratory data (Sui et al. 1994) and field observations, the parameters in the developed model are determined. Combined with the hydraulic model under the ice covered flow, the interaction between the discharge and the hydraulic pressure under the ice cover during the mechanical breakup of an ice cover is simulated. The impacts of the flow discharge during the mechanical breakup process on the dates of the breakup of a river ice cover are studied. Based on the field observations of breakups of the river ice cover in the Shiyaobu cross section of the Hequ Reach of the Yellow River, the hydraulic pressure on the ice cover under different discharges is determined by combing the model of the mechanical breakup of an ice cover with the model of river hydraulics under the ice covered flow. Results indicate that the simulated dates for the mechanical breakup of the river ice cover agreevery well with the field observations, although the calculated date for the river breakup in 1993 and 1994 is one-day earlier than the observed dates. It may be concluded that the proposed method for forecasting the date of the mechanical breakup of a river ice cover is applicable in practice.

Acknowledgement

This work was supported by the Hefei University of Technology (Grant No. GDBJ2008-020-Seed Grant for Ph. D.).

[1] SUI J., KARNEY B. and SUN Z. et al. Field investigation of frazil jam evolution–A case study[J]. Journal of Hydraulic Engineering, ASCE, 2002, 128(8): 781- 787.

[2] SUI J., KARNEY B. and FANG D. Variation in water level under ice-jammed condition–Field investigation and experimental study[J]. Nordic Hydrology, 2005, 36(1): 65-84.

[3] SUI J., WANG J. and HE Y. et al. Velocity profiles and incipient motion of frazil particles under ice cover[J]. International Journal of Sediment Research, 2010, 25(1): 39-51.

[4] WANG J., SUI J. and CHEN P. Numerical simulations of ice accumulation under ice cover along a river bend[J]. International Journal of Environmental Science and Technology, 2009, 6(1): 1-12.

[5] WANG J., SUI J. and GUO L. et al. Forecast of water level and ice jam thickness by using the back propagation neural network and support vector machine methods[J]. International Journal of Environmental Science and Technology, 2010, 7 (2), 215-224.

[6] SUI J., WANG J. and BALACHANDAR R. et al. Accumulation of frazil ice along a river bend[J]. Canadian Journal of Civil Engineering, 2008, 35(2): 158-169.

[7] SUI J., WANG D. and KARNEY B. Sediment concentration and deformation of riverbed in a frazil jammed river reach[J]. Canadian Journal of Civil Engineering, 2000, 27(6): 1120-1129.

[8] SUI J., HICKS F. and MENOUNOS B. Observations of riverbed scour under a developing hanging ice dam[J]. Canadian Journal of Civil Engineering, 2006, 33(2): 214-218.

[9] BELTAOS S. Threshold between mechanical and thermal breakup of river ice cover[J]. Cold Regions Scie- nce and Technology, 2003, 37(1): 1-13.

[10] BELTAOS S. Field studies of ice breakup and jamming in lower Peace River[J]. Cold Regions Science and Te- chnology, 2009, 56(2): 102-114.

[11] MAO Ze-yu, ZHAO Xue-feng and WANG Ai-min. Mechanism of longitudinal-crack formation during riverice breakup[J]. Advances in Water Science, 2009, 20(3): 434-437(in Chinese)

[12] YU D., YUSU Y. The research on the breakup of river ice cover prediction based on artificial neural network model[C]. 2011 IEEE 3rd International Conference on Communication Software and Networks. Xi’an, China, 2011, 17-19.

[13] ZHAO L., HICKS F. E. and FAYEK A. R. Applicability of multilayer feed-forward neural networks to model the onset of river breakup[J]. Cold Regions Science and Technology, 2012, 70(10): 32 - 42.

[14] ZHAO L., HICK F. E. and ROBINSON F. A. River breakup forecasting by hydro-meteorological data[C]. Proceedings, Annual Conference - Canadian Society for Civil Engineering. Ottawa, Canada, 2011, 3196- 3206.

[15] GUO Q. Applicability of criterion for onset of river ice breakup[J]. Journal of Hydraulic Engineering, ASCE, 2002, 128(11): 1023-1026.

[16] YANG Bo-yuan. Mechanics of materials[M]. Beijing, China: China Machine Press, 2002(in Chinese).

[17] LI Zi-jun, WANG Yong-xue and WANG Xi-wen et al. Effects of cement content and curing duration on physical and mechanical properties of synthetic ice for model test[J]. Journal of Hydraulic Engineering, 2004, (6): 41-45(in Chinese).

[18] WANG Fu-jun. Computational fluid dynamics[M]. Beijing, China: Tsinghua University Press, 2004(in Chinese).

10.1016/S1001-6058(11)60380-7

* Project supported by the National Natural Science Foundation of China (Grant No. 50979021), the Natural Science Foundation of Anhui Province (Grant No. 090415217).

Biography: WANG Jun (1962-), Male, Ph. D., Professor

SUI Jueyi, E-mail: sui@unbc.ca