SciELO - Scientific Electronic Library Online

vol.67 issue1Evaluation of weathering of basaltic lithotypes based on pH, electrical conductivity, and point-load strength testsFlexural behavior of hybrid steel fiber reinforced self-consolidating concretes author indexsubject indexarticles search
Home Pagealphabetic serial listing  

Services on Demand




Related links


Rem: Revista Escola de Minas

Print version ISSN 0370-4467

Rem: Rev. Esc. Minas vol.67 no.1 Ouro Preto Jan./Mar. 2014 



BEM numerical simulation of spillway flows


Simulação numérica de escoamentos sobre vertedores pelo MEC



Carlos Eduardo Ferraz de Mello

Graduate Program of Environmental Engineering - PROAMB. Federal University of Ouro Preto - UFOP/School of Mines - EM Department of Civil Engineering - DECIV.




The study of gravity-free surface flows presents difficulties, such as the nonlinearity of the dynamic boundary condition in the free surface, and also the fact that the location of this surface is not known a priori. Traditionally, this phenomenon has been investigated by physical models, but the progress of computer science and numeric methods has allowed more and more the successful use of mathematical models to simulate this type of flow. This work presents a boundary element method (BEM) numerical simulation of spillway flows with discontinuous linear elements. The solution procedure involves an iterative process in the determination of the free surface. The Newton-Raphson method is adopted together with the use of pseudo-nodes on the free surface and an empiric step factor (or damping factor) which controls the stability and the rate of convergence. An example of WES standard spillway shape is presented. The obtained results are compared with experimental data and they check the efficiency and good precision of the adopted method.

Keywords: Free surface flow, spillway, boundary element method.


O estudo de escoamentos com superfície livre, sujeitos à ação da gravidade, apresenta dificuldades como a não linearidade da condição de contorno dinâmica na superfície livre e, também, o fato da localização desta superfície não ser conhecida a priori. Tradicionalmente, esse fenômeno tem sido investigado por modelos físicos, mas os avanços da informática e dos métodos numéricos têm permitido cada vez mais o uso de modelos matemáticos com sucesso na simulação desse tipo de escoamento. Esse trabalho apresenta uma simulação numérica de escoamentos sobre vertedores pelo método dos elementos de contorno (MEC) com elementos lineares descontínuos. O procedimento de solução envolve um processo iterativo na determinação da superfície livre. É adotado o método de Newton-Raphson junto com o uso de pseudo-nós na superfície livre e de um parâmetro empírico de relaxação, que controla a estabilidade e a velocidade de convergência. Um exemplo de vertedor padrão tipo WES é apresentado. Os resultados obtidos são comparados com dados experimentais e comprovam a eficiência e a boa precisão do método adotado.

Palavras-chave: Escoamento com superfície livre, vertedor, método dos elementos de contorno.



1. Introduction

The application of numeric methods in water resource engineering is increasing more and more and in an accelerated way. In many cases, these techniques provide cheap and efficient alternatives for project verification and optimization, helping engineers to diagnose and solve possible problems.

Several water resource problems can be represented by the potential theory, such as sluice gates or spillways flows, examples of rapidly varied free surface flows governed by gravity. This kind of flow presents inherent difficulties for its solution. The difficulties, as pointed out by von Karman (1940) apud Cheng et al. (1981), arise not only from the nonlinear character of the dynamic boundary condition on the free surface, but also from the fact that the position of this surface is not known a priori. Since flows of this nature are characterized by highly curved contours, it is usually necessary, for a good numeric approach, the use of refined mesh or high order interpolation functions. Besides, due to the intrinsic nonlinearity of the problem, a great number of iterations seem to be inevitable in the solution procedure. Thus, the efficiency of the solution method is important and can become the main factor in the choice of the numeric method.

According to Cheng et al. (1981), serious numeric solutions to free surface flows began with Southwell & Vaisey (1946), using the finite difference method (FDM) with handmade calculations. Other finite difference solutions include those by McNown et al. (1955), Cassidy (1965) and, more recently, Liao et al. (1998) and Stelling & Zijlema (2003). Solutions with the finite element method (FEM) began with McCorquodale & Li (1971) for sluice gates. Among many other works that appeared in the literature using FEM, mentioned can be Ikegawa & Washizu (1973), Isaacs (1977), Betts (1979), Khan & Steffler (1996), Sankaranarayanan & Rao (1996) and Braess & Wriggers (2000). As examples of use of the boundary element method (BEM) in the solution of this problem type, mentioned can be Cheng et al. (1981), Liggett & Salmon (1981), Jovanovic (1987), Medina et al. (1991) and Machane & Canot (1997). Another adopted technique is the joint use of the complex function theory and integral equations solved numerically, as in the works of Merino (1996), Vanden-Broeck (1997) and Guo et al. (1998). Recently, the finite volume method (FVM) has also been used, as in the works of Olsen & Kjellesvig (1998), Unami et al. (1999), Song & Zhou (1999) and Qu et al. (2009). Still another technique recently adopted is the meshless particle scheme, and a good example is the Smooth Particle Hydrodynamics (SPH) formulation, as in the work of Ferrari (2010).

This work presents a numeric simulation of spillway flows through BEM classic formulation with discontinuous linear elements, having as basis the procedure proposed by Cheng et al. (1981).


2. Materials and methods

In the spillway flow problem here treated, only the irrotational steady state is considered. The stream function ψ formulation in a domain Ω is used and the governing differential equation is the Laplace's equation:

where x = (x,y) is a generic point of the domain. Both fixed boundaries and free surface are streamlines, therefore the value of ψ should be constant, which is:

being q the flow rate per unit width. On the free surface the dynamic boundary condition (Bernoulli's theorem) requires (Figure 1):

where: ν is the velocity.

g is the acceleration of gravity.

h is the free surface elevation measured from an arbitrary datum.

B is the constant of Bernoulli.

p is the pressure.

Starting from the stream function definition, the following relationship in the free surface is valid:

being n the unit normal from the free surface. Substituting Equation 5 in Equation 4 gives:

Either the flow rate, q, or the Bernoulli constant, B, are known, being the other obtained as part of the solution of the problem.

For the purpose of limiting a domain for the numeric solution, a cut is made at right angles to the primary velocity, at a certain distance upstream and downstream from the spillway crest. On the cut sections the following boundary condition is applied:

which means that there is no normal velocity to the main flow. This condition, although approximate, is applied sufficiently far away from the spillway crest so that it is quite accurate and any small error that occurs doesn't affect the interesting part of the flow (Cheng et al., 1981). The Figure 2 shows the boundary conditions on the problem domain.

The whole BEM discretization process and numeric approach to solve Equation 1 with boundary conditions of the Dirichlet type (Equations 2 and 3) and the Neumann type (Equation 7), which leads to a linear equation system in function of the unknowns y and n at the boundary functional nodes, is described in full detail in Mello (2003).

The inclusion of the free surface dynamic boundary condition (Equation 6) brings complications for the solution of the problem. This aspect is treated here in much the same way as it has been done for most of the finite difference and finite element solutions. That is, for problems with a known flow rate, q, a position of the free surface boundary is adopted and the problem solved from Equations 1, 2, 3 and 7. The Bernoulli constant is then calculated on the free surface by Equation. 6. If B is the same for all free surface points (nodes), the problem is solved. Otherwise, the adopted free surface is adjusted, iteratively, in order that the value of B becomes constant at all points. A similar procedure can be used for problems with a known Bernoulli constant, B. A position of the free surface is assumed and the problem is solved from Equations 1, 2, 6 and 7. The flow rate, q (Equation 3), is part of the solution. If q is constant for all free surface points, the problem is solved. Otherwise, an iterative process should be used to adjust the free surface elevation. The method of adjusting the free surface is problem dependent (Cheng et al., 1981).

The more common spillway flow problem is to solve for the Bernoulli constant, knowing the flow rate and determining the free surface profile. As indicated in Figure 1, there is a zone of uncertainty (Cheng et al., 1981) where a change in the water surface elevation has a negligible effect on the Bernoulli constant, B. Therefore, a direct iterative procedure supposing that the free surface points are independently adjusted cannot be applied. The corrections in the free surface elevation should be calculated simultaneously in all points. The method of Newton-Raphson is then adopted jointly with the concept of pseudo-nodes interpolated among the free surface points (real nodes). The following relationship is then obtained:

being {Δh} the elevation adjustment vector only for the real free surface nodes; {DB} is the departure of the Bernoulli constant for each free surface node (real and pseudo) from its objective value. For a given iteration k+1, the objective value of B is taken to be the average value of iteration k. The matrix [dB/dh] is computed by displacing, in turn, each real free surface node a small distance Dh from its adopted value (or the result of the previous iteration) and computing the variation in the Bernoulli constant at all free surface nodes, real and pseudo. This calculation requires a complete solution for each real free surface node, and the BEM does it in a quite efficient way. As a practical matter, a step factor (or damping factor), which controls the stability and the rate of convergence, is adopted in order that the adjustments to the free surface are always reasonable. Therefore, the values of the elevation in the iteration k+1 are obtained by the following equation:

where α is the step parameter. The value of this parameter is problem dependent and should be determined empirically. The solution is considered satisfactory when the calculated corrections on the free surface elevation become sufficiently small.


3. Results and discussion

In order to verify the efficiency of the developed procedure, some examples were simulated and the numeric results of one of them are presented. The chosen example is a standard WES type spillway designed for a flow rate q = 0.3744 m3/s.m and a height P = 1.52 m. This example allows validate the implemented formulation, since it has a known experimental result (Chow, 1959).

The domain was limited upstream to a distance (3.05 m) of twice the spillway height and downstream to a distance (1.16 m) in which the parament almost equaled the reservoir's bottom level. The domain discretization was made by 54 linear elements being 24 for the free surface. Discontinuous elements were used to considered the two values of the velocity ν in all the real corners.

The Tables 1 and 2 present, respectively, the main adopted and computed parameters for this example and Figure 3 shows the comparison of the free surface numeric result with the experimental water profile.

In relation to the data presented in Table 1, it is worthy to point out the following observations:

a) The step parameter a plays a fundamental role in the calculation march scheme, mainly in the initial iterations. To guarantee the stability in the initial phase of the process, it was almost always necessary to assume a small value. However, as a constant value along the whole procedure was adopted, it led to a slower convergence and a greater number of iterations. According to Jovanovic (1987), an automatic adjustment of this parameter at each iteration, through gradient type methods, would improve the rate of convergence.

b) The value of the displacement Dh to the numeric calculation of matrix [dB/dh] also has great influence in the success and efficiency of the iterative process. Very small values of Dh tend to increase the instability. Adopted was the same value of Dh for all free surface points and matrix [dB/dh] was calculated only once, in the beginning of the process, to make it more computationally efficient. According to Cheng et al. (1981), since matrix [dB/dh] changes little from one iteration to the next, it is not necessary to recalculate it each time. When it does have to be recomputed, it is only necessary to disturb and recompute those points in and near the uncertainty zone. As it is not clear when it is necessary to update the matrix, a more detailed study to optimize this procedure is suggested. Another possible change is the adoption of different values of the displacement Dh for each point, for instance, proportional to the difference between the value of the constant B in the point and an objective or average value of B.

c) The number of pseudo nodes has a certain importance, since in the case of linear elements, its absence generates instability. According to Liggett & Salmon (1981), the use of cubic spline elements eliminates the need of pseudo nodes. The procedure here developed allows the option of 2 and up to 9 pseudo nodes per element on the free surface, and was verified that starting from the value 4, no significant improvement of the results were found.

Now in relation to the data presented in Table 2, the following aspects can be highlighted:

a) The number of 59 iterations, for the circumstances commented previously, can be considered a good result. Cheng et al. (1981) didn't comment on this value in their work and Jovanovic (1987) reported a number of iterations about 160.

b) The value of the computed discharge coefficient, equal to 2.228, was very close to the experimental value 2.225 (Chow, 1959). This represents an error of only 0.13%. Cheng et al (1981) obtained an error of 0.50% and Jovanovic (1987) maximum mistakes about 3%.

c) The procedure also allows at its end the determination of the critical section by interpolation, since the value of the Froude number is calculated continually at all free surface points.

In relation to Figure 3, noticed can be a good adherence among the numeric results and the experimental profile (Chow, 1959), occurring just a small discrepancy in the area of greater curvature, where BEM slightly underestimated the surface water profile, meaning that the calculated velocities were a little larger. This can be associated to the fact that the mathematical model supposes that there are not energy losses along the flow over the spillway.


4. Conclusions

The example with known solution, here presented, confirms in general the applicability of the potential flow theory to the analysis of spillway flows, as well as validates the numeric simulation through BEM classic formulation with discontinuous linear elements, showing that BEM is an efficient numeric method to solve free surface potential flow problems.

The iterative procedure for the Newton-Raphson method adopted, and proposed by Cheng et al. (1981), led to very good results, but it still presents a slow convergence in the zone of uncertainty, taking a high number of iterations.

The role of the step parameter a, as well as the value of the displacement Dh, to the numeric calculation of matrix [dB/dh] should be further investigated in order to achieve an improvement on the rate of convergence of the iterative process.


5. Acknowledgments

The author wishes to thank FAPEMIG for the financial support.


6. References

BETTS, P. L. A variational principle in terms of stream function for free-surface flows and its application to the finite element method. Computers and Fluids, v. 7, p.145-153, 1979.         [ Links ]

BRAESS, H., WRIGGERS, P. Arbitrary lagrangian eulerian finite element analysis of free surface flow. Computer Methods in Applied Mechanics and Engineering, v.190, p. 95-109, 2000.         [ Links ]

CASSIDY, J. J. Irrotational flow over spillways of finite height. Journal of the Engineering Mechanics Division - ASCE, v. 91, n. EM6, p. 155-173, 1965.         [ Links ]

CHENG, A. H.-D., LIGGETT, J. A., LIU, P. L-F. Boundary calculations of sluice and spillway flows. Journal of the Hydraulics Division - ASCE, v. 107, n. HY10, p.1163-1178, 1981.         [ Links ]

CHOW, V. T. Open-Channel Hydraulics, New York: McGraw-Hill, 1959. 680p.         [ Links ]

FERRARI, A. SPH simulation of free surface flow over a sharp-crested weir. Advances in Water Resources, v. 33, p. 270-276, 2010.         [ Links ]

GUO, Y. et al. Numerical modelling of spillway flow with free drop and initially unknown discharge. Journal of Hydraulic Research, v. 36, n. 10, p. 785-801, 1998.         [ Links ]

IKEGAWA, M., WASHIZU, K. Finite element method applied to analysis of flow over a spillway crest. International Journal for Numerical Methods in Engineering, v.6, p. 179-189, 1973.         [ Links ]

ISAACS, L. T. Numerical solution for flow under sluice gates. Journal of the Hydraulics Division - ASCE, v. 103, n. HY5, p. 473-481, 1977.         [ Links ]

JOVANOVIC, M. B. Application of the boundary element method for the solution of spillway flows. In: BOUNDARY ELEMENT TECHNIQUES CONFERENCE (BETECH), 3. Proceedinga... Rio de Janeiro: 1987, Chap. Fluid Flow and Computational Aspects, p. 101-109.         [ Links ]

KHAN, A. A., STEFFLER, P. M. Modeling overfalls using vertically averaged and moment equations. Journal of Hydraulic Engineering, v. 122, n. 7, p. 397-402, 1996.         [ Links ]

LIAO, H., XU, W., WU, C. Residual pressure feedback method for simulation of free surface flow. Journal of Hydraulic Engineering, v. 124, n. 10, p. 1049-1052, 1998.         [ Links ]

LIGGETT, J. A., SALMON, J. R. Cubic spline boundary elements. International Journal for Numerical Methods in Engineering, v. 17, p. 543-556, 1981.         [ Links ]

MACHANE, R., CANOT, E. High-order schemes in boundary element methods for transient non-linear free surface problems. International Journal for Numerical Methods in Fluids, v. 24, p. 1049-1072, 1997.         [ Links ]

McCORQUODALE, J. A., LI, C.Y. Finite element analysis of sluice gate flow. Engineering Journal, v. 54, no. 3, p. 1-4, 1971.         [ Links ]

McNOWN, J. S., HSU, E-Y., YIH, C-S. Application of the relaxation technique in fluid mechanics. Transactions - ASCE, v. 120, p. 650-686, 1955.         [ Links ]

MEDINA, D. E. et al. A consistent boundary element method for free surface hydrodinamic calculations. International Journal for Numerical Methods in Fluids, v. 12, p. 835-857, 1991.         [ Links ]

MELLO, C. E. F. Boundary element method numerical analysis of sluice gate and spillway flows. Rio de Janeiro: PEC/COPPE/UFRJ - Federal University of Rio de Janeiro, 2003. (Doctoral Thesis) (In Portuguese).         [ Links ]

MERINO, F. Fluid flow through a vertical slot under gravity. Applied Mathematical Modelling, v.20, p. 934-939, 1996.         [ Links ]

OLSEN, N. R. B., KJELLESVIG, H. M. Three-dimensional numerical flow modelling for estimation of spillway capacity. Journal of Hydraulic Research, v. 36, n. 5, p.775-784, 1998.         [ Links ]

QU, J. et al. Numerical simulation of sharp-crested weir flows. Canadian Journal of Civil Engineering, v.3 6, p. 1530-1534, 2009.         [ Links ]

SANKARANARAYANAN, S., RAO, H. S. Finite element analysis of free surface flow through gates. International Journal for Numerical Methods in Fluids, v.22, p. 375-392, 1996.         [ Links ]

SONG, C. C. S., ZHOU, F. Y. Simulation of free surface flow over spillway. Journal of Hydraulic Engineering, v. 125, n. 9, p. 959-967, 1999.         [ Links ]

SOUTHWELL, R. V., VAISEY, G. Relaxation methods applied to engineering problems: XIII, Fluid motions characterized by free streamlines. Philosophical Transactions, Royal Society of London, Series A, v. 240, p. 117-161, 1946.         [ Links ]

STELLING, G., ZIJLEMA, M. An accurate and efficient finite-difference algorithm for non-hydrostatic free-surface flow with application to wave propagation. International Journal for Numerical Methods in Fluids, v. 43, p. 1-23, 2003.         [ Links ]

UNAMI, K. et al. Two-dimensional numerical model of spillway flow. Journal of Hydraulic Engineering, v. 125, n. 4, p. 369-375, 1999.         [ Links ]

VANDEN-BROECK, J.-M. Numerical calculations of the free-surface flow under a sluice gate. Journal of Fluid Mechanics, v. 330, p. 339-347, 1997.         [ Links ]



Artigo recebido em 02 de setembro de 2012.
Aprovado em 27 de novembro de 2013.

Creative Commons License All the contents of this journal, except where otherwise noted, is licensed under a Creative Commons Attribution License