Volume of Fluid (VOF) Method for Curved Free Surface Water Flow in Shallow Open Channel 

 

Li Ling , Chen Yongcan and Li Yuliang

Department of Hydraulic Engineering, Tsinghua University

Beijing 100084, China 

 

Abstract: Volume of fluid method is employed to predict the surface flow in a curved open channel. This method is based on the concept of a fractional volume of fluid (VOF), which is usually incorporated into the flow equations to track free water surface. The equations are discretized using the finite volume method. The results are agreed well with experimental data and VOF method can provide good results of water surface of the channel flow. It is shown that this method is more flexible and efficient than any other method for complicated free boundary configurations of curved channel.

1    INTRODUCTION

The computation of unsteady free-surface flow is required for the analysis of tidal oscillations, flows produced by the operation of control structures such as sluice gates, pumps or turbines, and flood waves generated by storms or failure of dams, dykes, or other structures. For such analyses, numerical models may be used as a management tool if the water level and wave speed can be accurately predicted over a wide range of conditions. Free boundaries are considered to be surfaces on which discontinuities exist in one or more variables. Examples are free surfaces, material interfaces, shock waves, or interfaces between fluid and deformable structures. More versatile and robust techniques for free surface flow modeling are the front-tracking methods . Front-tracking techniques are divided into two groups: surface-tracking and volume-tracking methods. In general the former gives a more accurate description of the free surface but the latter can handle complicated liquid regions more easily. Volume-tracking techniques define a tracer to follow the whole fluid region. The two commonly used techniques are the marker-and-cell(MAC) and the volume-of-fluid (VOF) techniques. In the marker-and-cell technique , hundreds of massless marker particles are added to the fluid. These particles are then advected in a Lagrangian sense using the average of Eulerian velocities in their vicinity. In the volume-of-fluid technique , a volume fraction parameter F is described for each of the Eulerian grid cells.A cell is assumed to be filled with liquid if F=1, empty if F=0 and partially full if 0<F<1. VOF based techniques can handle the most complex free surface flow problems. Surface breaking and merging can be treated with this technique, since the flow field calculations are decoupled from the free surface location identification. In present paper, we will provide the volume of fluid model (VOF) to solve the problems rising in the numerical treatment of free boundaries.

2    THE VOLUME OF FLUID MODEL

The VOF model can model two or more immiscible fluids by solving a single set of momentum equations and tracking the volume fraction of each of the fluids throughout the domain. Typical application includes the prediction of jet breakup, the motion of large bubbles in a liquid, the motion of liquid after a dam break, and the steady or transient tracking of any liquid-gas interface.

The VOF formulation relies on the fact that two or more fluids (or phase) are not interpenetrating. For each additional phase that you add to your model, a variable is introduced: the volume fraction of the phase in the computational cell. In each control volume, the volume fractions of all phases sum to unity. The fields for all variable and properties are shared by the phases and represent volume-averaged values, as long as the volume fraction of each of the phase is known at each location. Thus the variables and properties in any given cell are either purely representative of one of the phases, or representative of a mixture of the phases, depending upon the volume fraction values. In other words, if the qth fluid’s volume fraction in the cell is denoted as , then the following three conditions are possible:

             =0     the cell is empty (of the qth fluid)

             =1     the cell is full (of the qth fluid)

             0< <1   the cell contains the interface between the fluids

Based on the local value of , the appropriate properties and variables will be assigned to each control volume within the domain.

3    THE VOLUME FRACTION EQUATION

The tracking of the interface between the phases is accomplished by the solution of a continuity equation for the volume fraction of multi-phase flow. For the qth phase, this equation has the following form:

The volume fraction equation will not be solved for the primary phase; the primary-phase volume fraction will be computed based on the following constraint:

The properties appearing in the transport equations are determined by the presence of the component phase in each control volume. In a two-phase system, for example, if the phases are represented by the subscripts 1 and 2, and if the volume fraction of the second of these is being tracked, the density in each cell is given by

All other properties (e.g., viscosity) are also computed in this manner.

4    GOVERNING EQUATIONS

The Momentum Equation:

k-εturbulent model:

where the generated item of turbulent kinetic energy , the turbulent viscous coefficient  and the constant value: , , ,  , 0.09.

5    INTERPOLATION NEAR THE INTERFACE

In the paper, control-volume formulation requires that convection and diffusion fluxes through the control volume faces should be computed and balanced with source terms within the control volume itself. The geometric reconstruction scheme is used to calculate the face fluxes for the VOF model. The geometric reconstruction scheme represents the interface between fluids using a piecewise-linear approach. This scheme is the most accurate and is applicable for general unstructured meshes. It assumes that the interface between two fluids has a linear slope within each cell, and uses this linear shape for calculation of the advection of fluid through the cell faces.

6    NUMERICAL STABILITY CONSIDERATIONS

Numerical calculations have often computed quantities that develop large, high frequency oscillation in space, time, or both. This behavior is usually referred to as a numerical instability, especially if the physical problem being studied is known not to have unstable solutions. When the physical problem does have unstable solutions and if the calculated results exhibit significant variations over distances comparable to a cell width or over times comparable to the time increment, the accuracy of the results can not be relied on. To prevent this type of numerical instability or inaccuracy, certain restrictions must be adopted. In present paper, we use under-relaxation factors: under-relaxation factor is equal to 0.3 in the pressure equation; under-relaxation factor is equal to 0.5 in the momentum equation; under-relaxation factor is equal to 0.5 in the k andεequations.

7    CALCULATION OF UNSTEADY FLOW IN CURVED OPEN CHANNEL

The curved open channel is 924-meter long, 14-meter wide and 17-meter tall and comprises of linear slope segment, curvilinear slope segment, level segment and twisty segment To analyze unsteady flow in curved channel, Three-dimensional equations describing the conservation of mass, momentum and turbulent model are solved. Using VOF method tracks the free water-surface boundary.

Fig. 1    Curved open channel (unit:m.)

7.1    Boundary conditions

Numerical boundary conditions are applied at the inlet, exit and interface between water and air, and the curved channel walls. At the inlet, the velocity, water level and k, e are given; pressure is extrapolated from the interior solution. At the outflow boundaries the exit pressure is specified and the flow variables are evaluated using the zero-order Riemann invariant extrapolation. On the curved channel wall the no-slip boundary condition is applied. For unsteady flows, the pressure values on the curved channel wall are extrapolated from the interior solution.

7.2    Results and discussions

Table 1    Flow parameters of open channel

Case

Water level in reservoir (m)

Flux in channel (m3/s)

Inlet section of channel

Averaged water level (m)

Averaged velocity(m/s)

1

603.98

4026

11.70

24.58

2

598.82

3825

11.52

23.72

To verify the mathematical model presented in the previous sections, computed results are compared to laboratory test data. In the experiment, the measurement is done only in the middle x-y plane of the channel. Therefore, the experimental data are provided in the middle x-y plane that is situated at the middle of the channel.

Fig. 2    Pressure comparison of computed and experimental results for case 1

Fig. 3    Pressure comparison of computed and experimental results for case 2

Fig.2 shows the pressure comparison of computed and experimental results for Case 1, for which the flow depth at the entrance to the channel is 11.7m, Fig.3 shows the pressure comparison of computed and experimental results for Case 2, for which the flow depth at the entrance to channel is 11.52m. The experimental data are plotted using the dot and the computed results are plotted by continuous line.

Seen from the Fig.2 and Fig.3, good results of hydrodynamic pressure on the bottom of channel is given by the mathematical model . In the last portion of the channel, the flow was three-dimensional in the curved portion of the channel and there was a variation in the depth and hydrodynamic pressure across section of the channel. The mathematical model simulated the hydrodynamic pressure on the bottom of the channel satisfactorily. The flow in the channel was regular and had no neglect pressure on the bottom of the channel.

Fig. 4    Water level comparison of computed and experimental results for case 1

Fig. 5    Water level comparison of computed and experimental results for case 2

Fig.4 and Fig.5 shows the water level comparison of computed and experimental results for Case 1 and Case 2. The experimental data are plotted using the dot and the computed results are plotted by continuous line.

Seen from the Fig.4 and Fig.5, the computed results of water level are good agreement with the experimental data. For example, in the Case 2, water level of experimental data is 11.82m at the inlet of curvilinear segment and is 1.64m higher than one of computed result. The water level of experimental data is 6.18m at the outlet of curvilinear segment and is 0.34m higher than one of computed result. The water level of experimental data is 6.54m at the outlet of curved channel and is 0.12m higher than one of computed result. The water steadily flows in the linear slope segment and stably come into curvilinear segment, no exceptionable phenomenon occurring.

Fig. 6    Velocity comparison of computed and experimental results for case 1

Fig. 7    Velocity comparison of computed and experimental results for case 2

Fig.6 shows the velocity comparison for Case 1. Fig.7 shows the velocity comparison for Case 2. The velocity is measured on the point of average water depth in the middle x-y plane of curved channel. Seen from the two figures, the velocity is increasing slowly in the linear slope segment and increasing rapidly in the curvilinear segment. The computed results of water level are good agreement with the experimental data.

8    SUMMARY

The volume of fluid (VOF) technique has been presented as a simple and efficient means of numerically method for tracking free boundaries. It is used in the calculation of water flow in the three-dimensional curved channel. The computed results by the VOF method are well agreed with experimental data and VOF method can provide good results of water surface of the channel flow. The calculations for curved channel flow show that VOF method works extremely well and can handle the complex free surface. Here, knowing the initial fluid surface location and the velocity field, the surface is advected by using the volume fraction field, a surface is reconstructed based on the newly calculated volume fractions and finally the controlling equations are solved for the new fluid domain to find the velocity field. This procedure is repeated throughout the advection of the fluid region. There is no need for iteration between the flow field calculation and the interface location identification. It is shown that this method is more flexible and efficient than any other method for complicated free boundary configurations.

References

[1]    F.Mashayek, N.Ashgriz, A Hybrid finite-element-volume-of-fluid method for simulating free surface flows and interfaces, International Journal of Numerical Method in Fluids,1995,Vol.20,pp.1363-1380.

[2]    F.H.Harlow and J.F.Welch, Numerical calculation of time-dependent viscous incom pressible flow of fluid with free surface, Phys.Fluids,1965,Vol.8,pp.2182-2190.

[3]    H.Miyata and S.Nishimura, Finite-difference simulation of nonlinear ship waves,J.Fluid Mech.,1985,Vol.157,pp327-332.

[4]    H.Miyata, Finite-difference simulation of breaking waves, J. Comput. Phys.,1986, Vol.65, pp.179-187.

[5]    C.W.Hirt and B.D.Nichols, Volume of fluid (VOF) method for the dynamics of free boundaries, J. Comput. Phys.,1981,Vol.39,pp201-225.

[6]    Sun Degen and Cai Gongchun, the experimental research on open channel of Xiluodu Hydro-power Project, report of Hydraulics Laboratory, Hehai University, 1999.