Water Hammer in Branched Networks

 

MASSIMO GRECO, ARMANDO CARRAVETTA

 

Dipartimento di Ingegneria Idraulica ed Ambientale, Università di Napoli «Federico II» , Via Claudio 21, 80125 Napoli.

Tel.+39-081-7683427 Fax.+39-081-5938936 E-mail: grecom@unina.it.

 

 

ABSTRACT

Water hammer is investigated with reference to the simplest network composed of a single bifurcation. The equations describing the process are written in dimensionless form, and the parameters affecting their solution are singled out. Neglecting friction, and considering only instantaneous closures, seven parameters are necessary to identify the system considered. Numerical solutions, varying those parameters, have been computed, to check their influence on the maximum overpressures.

The results show that values of overpressure two times greater than the ones expected in a single pipeline can be found. Furthermore, drastic change in the computed overpressures can appear with small changes in diameter, wave speed or pipe length.

 

Keywords: Water hammer, pipe networks

 

INTRODUCTION

It is well known that water hammer in pipe networks may give rise to high pressures, due to the superposition of reflected pressure waves. This notwithstanding, the writers are not aware of previous detailed investigation of the process.

The problem is relevant nowadays for the common use of automatic valves in network regulation. When large variations in user demand happen often, as in distribution networks, water hammer occurs daily. In irrigation networks its effect was found to concur to the high rate of pipe failures those networks usually suffer [Brunone and Carravetta, 1994]. Actually, at least in Italy, water hammer has to be taken into account in the design of pipe networks, so it has been deemed useful to verify in depth its effects. Consider the simple network drawn in fig. 1.

Fig.1: Scheme of the system

 

It consists of three pipes and a node, marked as N. The first pipe, deriving from a constant level reservoir, has length L₁, cross section area A₁, and given properties so that its wave propagation velocity is a₁. In node N it splits its discharge Q₁ in two part, Q₂, and Q₃, flowing into the pipes 2 and 3. Subscripts 2 and 3 are attached to initial velocity V_o, length L, cross section area A, and wave propagation velocity a for the corresponding pipe. The reservoir level is h_s above the elevation datum.

Assume the following data: L₁=L₂=L₃=1000 m; v_o1=v_o3=0.9 m/s; v_o2=0; A₁=A₃=0.196 m²; A₂=0.098 m²; a₁=a₂=a₃=1000 m/s; h_s=150 m; and consider an instantaneous closure of the end valve of pipe 3. Neglecting friction, the hydraulic grade time history in node N is shown in fig. 2.

Fig. 2: Plot h(t) for L₁=L₂=L₃=1000 m.

 

In fig. 3 the computed time history is plotted assuming the same data as above, but varying L₂ to 800 m.

Fig. 3: Plot h(t) for L₁=L₂=1000 m, L₂=800m

 

A small change in one data item has produced a dramatic change in the computed pressures. Maximum pressures are larger, and, while they only last for very short time spans, no significant increase in pipe strength [Wakabayashi, 1986; Rusch, 1960] is to be expected for the short duration of the loading condition.

In this paper the water hammer equations for the very simple network drawn in fig. 1 are first written in non dimensional form, to reduce the number of intervening parameters. Then, the effects arising from many different combinations of the parameters have been computed. The results give the opportunity of some insights both on the numerical procedures to be used for the computations, and on the overall confidence that on such computations can be held.

 

THE MATHEMATICAL MODEL

For each pipeline, water hammer equations hold:

where v' is the mean velocity, h' the hydraulic grade, A the pipe cross section area, the fluid density, g the gravity acceleration, s' the space coordinate, t' the time. The friction term J can be assumed proportional to the square of velocity v', .

Solutions are usually obtained by numerical integration of the ordinary differential equations:

where a is the wave propagation velocity. Those equations can be derived by the method of characteristics under the assumption a>>V (Wylie et Al., 1993).

For the simple branched network of fig.1 considered in this study eq. (3) and (4) hold for each pipe, together with the continuity equation in node N, written, according to Wylie et al. (1992), as:

 

 

where is the summation of all instantaneous pipeline flows and a possible nodal flow, here taken always to be zero.

Being Q_oj, V_oj, and , respectively, rate of flow, velocity and friction in the initial steady state condition, and assuming:

 

 

equations (3) and (4), after substitution, become:

 

 

Considering, at first, only situations in which friction can be neglected (J₀₁=0) and applying the finite differences discretization with , where 1+n_j is the number of nodes in pipe j, the equations (5) and (6) can be written in the form:

 

where j=1..3 indicates pipe number, i is a generic node along the pipe j.

The continuity equation (7),using the equations (8) and (9), can be written as:

 

 

where:

 

Parameters P₂ and P₃ depend uniquely on parameters and already defined, and on the ratios and , through relations:

 

 

Restricting the analysis to an instantaneous closure of the valve at the end of pipe 3, no further equation is needed. Of course this is not a realistic simulation, but a crafted case study to gain insight into the problem. Then, under the restrictions introduced above, the problem is completely defined by the parameter set (I): . It seems useful to highlight that the parameter E defines the initial steady state conditions.

Although the choice of parameter set (I) is suitable to keep the equations (8-10) in a compact form, a different choice is to consider the alternative set of independent parameters (II): .The two parameter set, by equations (11) and (12) are completely equivalent, but a simpler physical interpretation is connected to parameter set (II); furthermore no mutual correlation exists among the parameters of set (II).

 

NUMERICAL SOLUTION

As said, usually in numerical modeling a first order finite difference method is applied. If the Courant number is 1, this method is exact. Otherwise, interpolation along time or space is needed. Interpolation introduces a significant numerical diffusivity, that smoothens the sharp pressure fronts that actually happen [Streeter, 1972; Ghidaoui and Karney, 1993; Karney and Ghidaoui, 1997]. With instantaneous interception of the flow, a superposition of square waves represents the solution in every point in the space-time domain. Whichever finite degree of interpolation will filter out some of the higher frequencies, and will so dampen the pressure peaks, diffusing them.

In a single pipe it is easy to select Ds and Dt so that the Courant number is 1. It may be impossible to do so when more than one pipe is involved in computations, if lengths and wave propagation velocities are arbitrary, and a single Dt value is sought for the entire system. Otherwise, by allowing the selection of a different Dt for each pipe, interpolation will be needed when dealing with the node equation [Wiggert and Sundquist, 1977]. From the practical simulation point of view, in order to avoid interpolations, it has been suggested sometimes to alter slightly the pipe lengths, or, better, the wave propagation velocities [Chaudry, 1987; Wylie et Al., 1993]. Many Authors maintain that interpolation cannot cause too much inaccuracy to calculations if maneuvers are slow enough as is usually the case in real hydraulic systems. Furthermore, higher order interpolation schemes have been proposed, like cubic splines, hermite polynomials and so on [Shimada and Okushima, 1984; Sibetheros et Al, 1991].

Sometimes, also when it is possible to select a Dt time step so that the Courant numbers are 1 in all the pipes, it is inconvenient to do so for the size of the needed time step is so small that an unusual amount of computations will be required.

The effects of interpolation have been tested considering the example given in the introduction. Computations have been performed with a Courant number equal to 0.9 and results are shown in fig 4. As can be seen by comparison with fig. 3, interpolating schemes severely underestimate the maximum pressures. Lower values of the Courant number make things worse.

Fig. 4: Plot h(t) for L₁=L₂=L₃=1000 m (with interpolation)

 

Therefore, all calculations have been performed always assuming a Courant number equal to 1 and a single time step Dt for the whole system. Some of the parameters of the network, within the parameter set (II), have been varied and their influence on the computed pressures will be discussed in detail. The range of variation of the parameters are reported in Table I. For completeness in the same table are reported also the corresponding ranges of the parameters of set (I).

 

 

PARAMETERS INFLUENCE

Maximum computed pressures in the node N for all simulations are shown in Table II.

Table II
	Max h	N (%) h<1	N (%) h<1.3	N (%) h<1.6	N (%) h<1.9
E=0.0	2.00	30.44	81.13	97.75	99.91
E=0.5	1.00	99.95	100.0	100.0	100.0

 

In the same table are reported the percentage of cases with h less than 1-1.3-1.6-1.9. For E=0.0 the maximum value of h was 2.00; in the 50 % of cases the computed pressures were between 1.0 and 1.3. For E=0.5 the values of h were almost always less than 1.00.

Then, limiting our analysis to the case E=0.0, the values of h are plotted as a function of the parameter in fig.5.

Fig. 5: Plot of h vs. E=0.0).

 

No regularity is easily spotted in the numerical results. In Table III are reported the maximum values of h found for the different values of .

 

 

Table III
Case						h
A	0.6	1.0	0.6	1.000	0.833	1.12
B	0.6	1.0	0.5	1.000	0.833	1.56
C	0.8	0.5	0.5	1.087	1.000	1.30
D	0.7	0.5	0.5	1.087	1.000	1.91
E	1.0	1.0	0.5	1.000	1.000	0.81
F	1.0	1.0	0.5	1.042	1.000	1.54

 

These values correspond to the larger values plotted in fig.5 for each . The highest value, h=2.0, was found for the value . For the smaller values of , i.e. for larger wave speed of pipe 3 compared to pipe 1, the values of h were larger than 1.0, independently from the value of the other parameters. For the larger values of the values of h were variable approximately between 0.6 and 1.6, but no continuous dependence of results on parameter seems to exist, as is typical of chaotic systems.

The values of h, found in all the calculations with E=0.0 and , are plotted in fig.6 as a function of .

Fig. 6: Plot of h vs. (E=0.0: =1)

 

For each value of , the smaller values of h are the ones obtained in the symmetric case, with (square symbols). The maximum values of h, found in each set of calculations with a fixed value of , seems to decrease for increasing values of , i.e. decreasing the diameter of pipe 3 compared to pipe 1.

More interesting is to compare (Table IV) the results obtained from computations in three cases, where only a small difference in one of the parameters exists. Considering cases A and B the value of changes by 20.0% while the value of h changes from 1.12 to 1.56 (39%). In cases C and D the value of changes by 14.3% while the value of h changes from 1.30 to 1.91 (47%). In cases E and F the value of changes by only a 4.2% while the value of h rises from 0.81 to 1.54 (90%).

The difference in the parameters considered above may be less than the actual uncertainty, at design time, on the exact parameter value for a real plant. The variation in the maximum pressures is large enough to worry about. Therefore, we suggest that, in performing water hammer computations as part of the design procedure for real networks, more than one value should be considered for parameters that are subject to uncertainty.

Of course including friction in the computations reduces the amplitude of superpositions that occur late in the process. As an example the same two cases of figures 2 and 3 have been recomputed assuming a quadratic friction law, with a steady state head loss along the pipeline equal to 5, 10, 25 % of the static head. The maximum values, h, for the case of fig. 2 and of fig.3 are reported in Table IV.

 

 

Even more significant is to consider closure maneuvers that last more than 1 rhythm. Again for the cases of figure 2 and 3, without considering friction, the results are reported in table IV for closure T=1, 2, 3, rhythms.

 

CONCLUSIONS

Bifurcations in pipe systems are certainly relevant to the evolution of water hammer transients. A systematic investigation on the influence that branching geometry has on the value of overpressures during the transient has been presented here.

The parameters influencing the phenomenon have been evidenced and calculations performed while varying the parameters. The maximum value of overpressure , computed during each simulation in the node of bifurcation, N, was recorded.

The largest value of h was two, corresponding to an overpressure twice the value expected from the sudden closure in a single pipe. The comparison between the h values in different simulation show no continuous dependency on the parameters.

 

Finally, the writers want to stress that differences in the values of pipe diameter, wave speed or length not exceeding the approximations usually taken at design time can lead to drastic changes in the overpressure. Therefore, much care should be exercised when these data are known only approximately or when predictions are computed changing length or wave speed in order to avoid interpolation.

 

REFERENCES

Brunone B., Carravetta A. (1994): Daily management of irrigation pipe networks - Water hammer pressure surges and happened damages, 2^nd International Conference on Pipeline System, BHRA Group.

Chaudry M.H. (1987): Applied Hydraulic Transients, 2^nd Edition, Van Nostrand Reinhold Company Ed., New York, USA.

Ghidaoui M.S., Karney B.W. (1994): Equivalent Differential Equation in Fixed-Grid Characteristics Method, Journal Hydraulic Engineering, ASCE, Vol.120, N.10.

Karney B.W., Ghidaoui M.S. (1997): Flexible Discretization Algorithm for Fixed-Grid MOC in Pipelines, Journal Hydraulic Engineering, ASCE, Vol.123, N.11.

Liou C.P. (1983): Calculation of transients in batched pipelines, 4^th International Conference on Pressure Surges, BHRA Group.

Rusch H.(1960): Researches toward a general flexural theory for structural concrete, Journal of the American Concrete Institute, Vol.32, N.1

Shimada M., Okushima S. (1984): New Numerical Model and Technique for Waterhammer, Journal Hydraulic Engineering, ASCE, Vol.110, N.6.

Streeter V.L. (1972): Numerical Methods for Calculations of Transient Flow, 3^th International Conference on Pressure Surges, BHRA Group

Sibetheros I.A., Holley E.R., branski J.M. (1991): Spline Interpolation for Water Hammer Analysis, Journal Hydraulic Engineering, ASCE, Vol.117, N.10.

Wakabayashi M. (1986): Design of Earthquake Resistant Buildings , McGrawHill.

Wiggert D.C., Sundquist M.J. (1977): Fixed-Grid Characteristic for Pipeline Transients, Journal of the Hydraulic Division, ASCE, Vol.103, N.HY12.

Wylie E. B., Streeter V.L., Suo L. (1993): Fluid Transient in Systems, Prentice-Hall Ed., USA