This paper develops quadratic approximations to the classical models of Darcy–Weisbach and Hazen–Williams for frictional head loss in water pipes. Key elements of the technique are selecting the approximation domain and minimizing the relative error. A comparison for individual pipes over a range of diameter, roughness, and Reynolds number showed that the approximation provides accuracy consistent with the experimental error in the classical equations. Two benchmark water distribution networks are also considered. In these systems, flows and pressures computed using the approximation were consistent with simulations using the classical models. Applications of the approximation include mathematical optimization problems where polynomial expressions are desirable.

NOTATION

  • a

    dimensional coefficient

  • b

    dimensional coefficient

  • f

    Darcy friction factor [–]

  • g

    gravitational acceleration [m/s/s]

  • hf

    frictional head loss [m]

  • approximation of frictional head loss [m]

  • ks

    roughness height [m]

  • AE

    constants used to compute a and b for Hazen–Williams

  • C

    Hazen–Williams C value

  • D

    pipe diameter [m]

  • F

    integral of relative errors

  • L

    pipe length [m]

  • Q

    flow rate [m3/s]

  • Q1

    flow rate at lower end of approximation interval [m3/s]

  • Q2

    flow rate at upper end of approximation interval [m3/s]

  • Re

    Reynolds number

  • V

    average fluid velocity [m/s]

  • α

    constant part of Hazen–Williams formula with respect to flow

  • ε

    relative error

INTRODUCTION

The equations describing frictional head loss in water pipes are well known and widely applied. Head loss is described by the Darcy–Weisbach equation, or in certain flow regimes using empirical formulas such as those of Hazen–Williams (Liou 1998; Christensen et al. 2000). Today's network modeling software packages such as EPANET (Rossman 2000) implement these approaches and are widely used for planning and operation of water systems. However, applications for network models are becoming more sophisticated to include data assimilation and real-time modeling (Hutton et al. 2014), optimization (Broad et al. 2010), and inverse problems. In these applications, the properties of the classical head loss equations with respect to complexity, differentiability, and smoothness limit the techniques that can be applied.

Techniques for solving head loss equations depend on the problem studied and friction loss model. Applications in network simulation often use a Newton-type algorithm that requires a Jacobian matrix. The Jacobian can be expensive to compute particularly for the Darcy–Weisbach model (Simpson & Elhay 2011). With the Hazen–Williams equations, the Jacobian can become singular near-zero flow causing simulations to fail (Elhay & Simpson 2011). In mathematical optimization problems, conservation of energy along pipes often appears as a constraint. The challenge of using the classical equations as constraints is that the rational exponents appearing in the Hazen–Williams formula or in explicit approximations for Darcy–Weisbach may have singularities at zero flow. A quadratic approximation for head loss allows to formulate a variety of network optimization problems as bounded polynomial optimization problems, non-convex problems whose global optima can be approximated by sequences of semi-definite optimization problems, that is, a class of convex optimization problems which can be solved efficiently (Lasserre 2001).

Numerous relationships for pipe head loss have been proposed. The empirical formulas of Hazen and Williams and Manning give explicit equations, but apply over a limited range of flows. In contrast, the Darcy–Weisbach equation for head loss is accurate and generally applicable, but requires a friction factor that depends on the flow regime. For turbulent flows the factor may be computed from the implicit Colebrook–White formula, or in the limit of fully rough flow from the explicit formula of Prandtl and Karman (White 1999). Many explicit approximations of the Darcy friction factor are also available, most notably those of Swamee & Jain (1976) for popularity and Sonnad & Goudar (2006) for accuracy. A recent summary of approximations for the friction factor is given by Giustolisi et al. (2011).

The need to simplify the Darcy–Weisbach formula for use in optimization problems has been recognized by a few earlier workers. Valiantzas (2008) proposes an explicit power law approximation that provides accuracy of ±4.5%. The fractional exponents in this approximation cause the same numerical difficulties as the Hazen–Williams formula. Burgschweiger et al. (2009) propose a globally smooth explicit approximation for pipe head loss that includes square root functions. The approximation is asymptotically consistent with classical equations, but as will be shown subsequently, may underestimate friction losses in certain situations. Gleixner et al. (2012) use the Prandtl–Karman relation which can also underestimate friction losses and has a singularity in the inverse derivative at zero flow.

The contributions of the present paper are quadratic approximations for friction losses in water pipes that minimize the relative error and apply to water systems parameterized by either the Darcy–Weisbach or Hazen–Williams models. The primary motivation for this approximation is to model head loss with high accuracy while providing a mathematical form, i.e., a polynomial approximation with integer powers, which is well-behaved when embedding it in constraints of optimization problems over water networks. Some examples of such problems include network design (Alperovits & Shamir 1977), valve placement (Reis et al. 1997), pump scheduling (Sterling & Coulbeck 1975), and valve setting (Vairavamoorthy & Lumbers 1998) among others. In particular, the work of Bragalli et al. (2011) on the network design problem addresses the numerical difficulties of rational exponents in the head loss function by fitting a quintic near-zero. An additional effect of the proposed approximation is a lower complexity for computing the head loss function and its derivative. The reduced complexity may also prove beneficial where hydraulic simulators make many runs guided by genetic algorithms (Savic & Walters 1997), simulated annealing (McCormick & Powell 2004), and other heuristics used for water network optimization.

BACKGROUND

The head loss due to friction in water distribution systems is well studied and may be described by the Darcy–Weisbach equation (White 1999) 
formula
1
where f is the Darcy friction factor; L is the pipe length; D is the diameter; V is the average fluid velocity; and g is the gravitational constant, all in consistent units.

The method for computing the Darcy friction factor varies according to whether the flow is laminar or turbulent. In the laminar regime, Poiseuille flow occurs and the friction factor is , where is the Reynolds number, and is the kinematic viscosity. Thus, under laminar conditions .

In the turbulent regime, the friction factor's behavior with respect to varies according to flow conditions as diagramed by Moody (1944). For fully rough turbulent flow, the factor is given by the classical formula of Prandtl–Karman 
formula
2
where ks is the roughness height and D is the diameter. In Equation (2), f is constant with respect to Re and thus . In turbulent flows that are not fully rough f depends in a complex way on Re. In either case, the value of f may be found graphically from the Moody diagram, or by inverting the implicit relationship of Colebrook (White 1999) 
formula
3
The value of f may also be estimated using one of the many explicit approximations in the literature (see Giustolisi et al. (2011) for a recent summary).
In the transition from laminar to fully rough flow, the Hazen–Williams formula also applies. In SI units, friction loss based on the Hazen–Williams formula is computed from American Water Works Association (2004), 
formula
4
where hf is the head loss due to friction in meters; C is the Hazen–Williams coefficient; D is the inside diameter of the pipe in meters; L is the length of pipe in meters; and Q is the discharge in cubic meters per second. The Hazen–Williams formula applies only to water flowing in pipes 5 cm in diameter or larger and at speeds less than 3 m/s (Franzini & Finnemore 1997). A more detailed discussion of the region of application for the Hazen–Williams formula is given by Liou (1998) and discussions.

Although Equations (3) and (4) are often treated as exact, they are fitted to experimental results and therefore contain some error. Within the flow range where tests were performed, Liou (1998) reports errors from −15 to 10% in the Hazen–Williams formula depending on pipe material. On the Darcy–Weisbach side, the Colebrook equation for f has an average error of 11% with respect to Colebrook's data (Yoo & Singh 2005).

QUADRATIC APPROXIMATION FOR PIPE FRICTION

In some applications such as real-time simulation or mathematical optimization, a simpler relationship than Equation (1) with Equation (3), or an explicit approximation for f, that provides a sufficiently accurate representation of the underlying physics is desirable. The limiting case of fully rough flow suggests a squared function. Indeed, the head-flow curve for a pipe shows strongly quadratic behavior in the turbulent range. Thus, the following quadratic approximation for the pipe head loss curve is proposed: 
formula
5
where a and b are unknown dimensional coefficients and Q is the volumetric flow rate. The intercept is dropped to enforce the origin as a point on the curve. The linear term is retained to provide a degree of freedom for improving the fit, and in recognition of the fact that turbulent flows outside of the fully rough range are not purely a squared function of the flow.
The derivative or reciprocal derivative of the head loss function is needed in simulation and optimization applications. Using a quadratic, the derivative expression is easily found and well behaved 
formula
6
Notice that from the left and the right. The reciprocal of the derivative is non-zero provided for .

In developing an approximation function one must consider several factors:

  • the functional form of the approximation;

  • a merit function for evaluating candidates;

  • the relevant range of the original function;

  • the method for finding coefficients.

This work uses the functional form of Equation (5) whether Darcy–Weisbach or Hazen–Williams is used as an underlying model. Coefficients are selected to minimize the ‘relative error’ (MRE) over the approximation range. When the coefficients of Equation (5) are found in this way, the approximation is termed a ‘MRE quadratic’. The exact merit function, approximation range, and method of finding coefficients differ slightly between head loss models as discussed below.

Darcy–Weisbach

The preferred method for modeling head loss on water networks is the Darcy–Weisbach law of Equation (1). If the friction factor f were constant the head loss curve would be purely quadratic; this is the primary motivation for developing a quadratic approximation. The fit proposed here minimizes the sum of squared relative errors over the turbulent flow range . Owing to the implicit form of Colebrook's equation for f, and the complexity of explicit approximations, the fit is obtained numerically by solving the least squares problem 
formula
7
By using as the weight for each point, the coefficients minimize relative rather than absolute errors. A vector of flow rates having length N defines the approximation interval. The lower end of the approximation interval corresponds to . At the upper end of the interval, QN corresponds to the maximum flow rate expected in the application based on a preliminary examination of the network. A sufficiently large number of points should be used so that the values of the fitted coefficients a and b do not depend on the number of points.

A quadratic fit to the Darcy–Weisbach head loss curve in a sample pipe is shown in Figure 1. The sample pipe has length 1,000 m, diameter 100 mm, and roughness height 0.3 mm. Considering water with a kinematic viscosity of 1.002 × 10−6 m2/s, a MRE quadratic approximation using 2,000 points over the range yielded coefficients and for use with flow rates in L/s. The relative error ranges from −1.9 to 8% with the largest relative errors occurring at low flow rates.

Figure 1

MRE quadratic approximation for head loss according to the Darcy–Weisbach formula for an example pipe where L = 1,000 m, D = 100, and ks = 0.3 mm. The fit used an interval of .

Figure 1

MRE quadratic approximation for head loss according to the Darcy–Weisbach formula for an example pipe where L = 1,000 m, D = 100, and ks = 0.3 mm. The fit used an interval of .

Because the quadratic fit is developed for a particular flow range, the accuracy of the approximation outside of the range is of interest. In the laminar region, head losses are linear and the quadratic approximation over the turbulent range does not provide a very good estimate (Figure 2). However, the errors noted here do not affect the flow distribution in the networks studied because head losses are low. Figure 2 also illustrates why the transition flow range is excluded from the approximation: trying to minimize relative errors through the jump in the transition range sacrifices accuracy everywhere else along the curve. On the upper end of the flow range, the shape of the relative error curve suggests good agreement with the theoretical curve beyond the range selected for approximation. In the example of Figure 1, applying the approximation at Re = 2 × 105, twice the fitted value, incurs an error of 2.2%.

Figure 2

Head losses at low flows for an example pipe where L = 1,000 m, D = 100, and ks = 0.3 mm. The quadratic approximation uses an interval of .

Figure 2

Head losses at low flows for an example pipe where L = 1,000 m, D = 100, and ks = 0.3 mm. The quadratic approximation uses an interval of .

Visual inspection of the Moody diagram shows that the friction factor varies over a smaller range as the relative roughness of the pipe increases. This smaller range implies that rougher pipes, such as those in aging distribution systems, are well approximated by a quadratic function. Indeed the fit improves as pipes become rougher (Figure 3(a), 3(c), and 3(e)).

Figure 3

MRE quadratic approximation for head loss at average velocities up to 3 m/s. In these pipes L = 100 m, D = 100, 200, and 300 mm; and roughness values are noted on the plots. (a), (c), and (e) use Darcy–Weisbach and (b), (d), and (f) use Hazen–Williams.

Figure 3

MRE quadratic approximation for head loss at average velocities up to 3 m/s. In these pipes L = 100 m, D = 100, 200, and 300 mm; and roughness values are noted on the plots. (a), (c), and (e) use Darcy–Weisbach and (b), (d), and (f) use Hazen–Williams.

Also shown in Figure 3 is the curve given by Prandtl–Karmman's formula and the approximation due to Burgschweiger et al. (2009). The Burgschweiger approximation is given by 
formula
8
where the constants r, b, and c are found from pipe properties and the values of a and d are chosen to match a desired slope at . In the comparisons shown here, a and d are chosen to match the slope at which occurs under laminar flow according to the law of Hagen–Poiseuille. For pipes with roughness heights of 10 and 1 mm, both a quadratic approximation, the law of Prandtl–Karman, and Burgschweiger's approximation fall directly on the ideal curve. As a roughness height decreases, deviation of the approximations from the ideal curve becomes visible.

Hazen–Williams

The Hazen–Williams formula for pipe head loss applies over a limited range; outside this range it does not reflect the underlying physics and has considerable errors compared to experimental results (Liou 1998; Franzini & Finnemore 1997). These deficiencies notwithstanding, the equation remains in wide use because the explicit formulation is easy to compute. Since existing models of real networks often use Hazen–Williams, a quadratic approximation is provided here.

A quadratic approximation to the Hazen–Williams formula can be found numerically as described above for Darcy–Weisbach. However, with its explicit functional form a quadratic approximation to the Hazen–Williams formula can also be obtained analytically by minimizing the integral of relative errors 
formula
9
where and specify the range of interest and α is the part of Equation (4) that does not vary with flow rate . The integral of Equation (9) is found and the resulting form is differentiated by a and b. These differentials are set equal to zero, giving a 2 × 2 system of linear equations. Solving the system yields expressions for a and b as functions of the approximation interval 
formula
10
where 
formula
11
A quadratic approximation to the Hazen–Williams curve for a sample pipe is shown in Figure 4 along with the relative error curve. Several points of interest are noted:
  1. the quadratic fit provides errors within 10% over the fitted range of flows;

  2. since the relative error is minimized, there is a singularity in the relative error at zero flow;

  3. the shape of the relative error curve suggests a method for choosing the approximation interval.

Figure 4

MRE quadratic approximation for head loss according to the Hazen–Williams formula for an example pipe where L = 1,000 m, D = 100 mm, and C = 120. The fit used an interval of 0.1–8 L/s where the lower bound was chosen to give a minimum relative error of 10%.

Figure 4

MRE quadratic approximation for head loss according to the Hazen–Williams formula for an example pipe where L = 1,000 m, D = 100 mm, and C = 120. The fit used an interval of 0.1–8 L/s where the lower bound was chosen to give a minimum relative error of 10%.

Choosing the approximation interval is a key step because the values of a and b depend on the interval. The value of should be chosen to reflect the maximum flow likely to be encountered in the application. A value can be chosen through investigation of the network or selection of a maximum velocity or power dissipation in each pipe. Once a value for is established, may be selected so that the minimum of the relative error curve sits at an error tolerance. The location of the minimum relative error is found from calculus 
formula
12
With Equation (12) and an error tolerance, εtol, is found at the root of 
formula
13
This approach to setting the approximation interval is conservative in that the worst underestimate of head losses is specified. Head losses which are overestimated occur only at the extreme ends of the interval and relative errors are largest for small flows. A quadratic approximation to the Hazen–Williams head loss curve is shown for a range of diameters, flows, and roughnesses in Figure 3(b), 3(d), and 3(f).

COMPARISON ON PIPE NETWORKS

A numerical model was created to study the accuracy of the MRE quadratic approximation for networks of pipes. The model solves the set of nonlinear equations for steady-state network hydraulics using the global gradient algorithm of Todini & Pilati (1988). For the MRE quadratic, Jacobian terms are given by Equation (6). For the Darcy–Weisbach model, the Jacobian matrix given by Simpson & Elhay (2011) is used. Note that the coefficients of the MRE quadratic approximation are pre-computed, and so are treated as constants by the simulator. Because the goal here was a comparison of accuracies, the implementation was not optimized for speed. In the following subsections, the simulation model is applied on two networks from the literature.

Twenty-five node network

The first network considered here is the benchmark network of Sterling & Bargiela (1984) that has also been studied by Jowitt & Xu (1990), Reis et al. (1997), Vairavamoorthy & Lumbers (1998), and Nicolini & Zovatto (2009) among others. The network layout is shown in Figure 5. With three sources, 22 junctions and 37 links, the network is highly connected. Head loss parameterization is given using Hazen–Williams coefficients.

Figure 5

Twenty-five node network due to Sterling & Bargiela (1984).

Figure 5

Twenty-five node network due to Sterling & Bargiela (1984).

Quadratic approximations for head loss in each pipe were developed using the methods given above with parameter corresponding to a velocity of 2 m/s and . For convenience of other workers, network information including coefficients of a quadratic approximation is provided in Tables 1 and 2.

Table 1

Node elevations and demands for network due to Sterling & Bargiela (1984) 

ID Elev [m] Demand [L/s] 
18 
18 10 
14 
12 
14 30 
15 10 
14.5 
14 20 
14 
10 15 
11 12 10 
12 15 
13 23 
14 20 
15 20 
16 10 
17 
18 
19 10 
20 
21 10 
22 15 20 
23 54.66 – 
24 56.4 – 
25 54.5 – 
ID Elev [m] Demand [L/s] 
18 
18 10 
14 
12 
14 30 
15 10 
14.5 
14 20 
14 
10 15 
11 12 10 
12 15 
13 23 
14 20 
15 20 
16 10 
17 
18 
19 10 
20 
21 10 
22 15 20 
23 54.66 – 
24 56.4 – 
25 54.5 – 
Table 2

Link properties including coefficients of minimum relative error quadratic approximation for network due to Sterling & Bargiela (1984) 

ID Node 1 Node 2 Length [m] Diameter [m] Roughness a [s2/m5b [s/m2
23 606 0.457 110 60.689 0.40258 
23 24 454 0.457 110 45.467 0.30160 
24 14 2,782 0.229 105 10,818 18.018 
25 14 304 0.381 135 53.365 0.24605 
10 24 3,383 0.305 100 3271.2 9.6655 
13 24 1,767 0.475 110 144.92 1.0385 
14 13 1,014 0.381 135 178.00 0.82070 
16 25 1,097 0.381 61,494 283.53 
1,930 0.457 110 193.28 1.282 
10 5,150 0.305 10 354,180 1046.5 
11 12 13 762 0.457 110 76.312 0.50622 
12 15 16 914 0.229 125 2573.3 4.2861 
13 17 16 822 0.305 140 426.24 1.2594 
14 18 17 411 0.152 100 14563 10.686 
15 20 18 701 0.229 110 2500.8 4.1654 
16 19 17 1,072 0.229 135 2617.2 4.3593 
17 20 19 864 0.152 90 37,211 27.307 
18 21 20 711 0.152 90 30,621 22.471 
19 21 15 832 0.152 90 35,832 26.295 
20 22 15 2,334 0.152 100 82,701 60.689 
21 12 15 1,996 0.229 95 9341.8 15.560 
22 11 12 777 0.229 90 4019.6 6.6952 
23 10 11 542 0.229 90 2803.9 4.6702 
24 12 1,600 0.457 110 160.23 1.0629 
25 10 249 0.305 105 219.97 0.65000 
26 443 0.229 90 2291.7 3.8172 
27 743 0.381 110 190.59 0.87872 
28 22 931 0.229 125 2621.1 4.3659 
29 22 21 2,689 0.152 100 95,280 69.920 
30 326 0.152 100 11,551 8.4767 
31 844 0.229 110 3011.0 5.0151 
32 1,274 0.152 100 45,142 33.127 
33 1,115 0.229 90 5768.1 9.6076 
34 615 0.381 110 157.75 0.7273 
35 22 1,408 0.152 100 49,890 36.611 
36 500 0.381 110 128.25 0.59133 
37 300 0.229 90 1552.0 2.5850 
ID Node 1 Node 2 Length [m] Diameter [m] Roughness a [s2/m5b [s/m2
23 606 0.457 110 60.689 0.40258 
23 24 454 0.457 110 45.467 0.30160 
24 14 2,782 0.229 105 10,818 18.018 
25 14 304 0.381 135 53.365 0.24605 
10 24 3,383 0.305 100 3271.2 9.6655 
13 24 1,767 0.475 110 144.92 1.0385 
14 13 1,014 0.381 135 178.00 0.82070 
16 25 1,097 0.381 61,494 283.53 
1,930 0.457 110 193.28 1.282 
10 5,150 0.305 10 354,180 1046.5 
11 12 13 762 0.457 110 76.312 0.50622 
12 15 16 914 0.229 125 2573.3 4.2861 
13 17 16 822 0.305 140 426.24 1.2594 
14 18 17 411 0.152 100 14563 10.686 
15 20 18 701 0.229 110 2500.8 4.1654 
16 19 17 1,072 0.229 135 2617.2 4.3593 
17 20 19 864 0.152 90 37,211 27.307 
18 21 20 711 0.152 90 30,621 22.471 
19 21 15 832 0.152 90 35,832 26.295 
20 22 15 2,334 0.152 100 82,701 60.689 
21 12 15 1,996 0.229 95 9341.8 15.560 
22 11 12 777 0.229 90 4019.6 6.6952 
23 10 11 542 0.229 90 2803.9 4.6702 
24 12 1,600 0.457 110 160.23 1.0629 
25 10 249 0.305 105 219.97 0.65000 
26 443 0.229 90 2291.7 3.8172 
27 743 0.381 110 190.59 0.87872 
28 22 931 0.229 125 2621.1 4.3659 
29 22 21 2,689 0.152 100 95,280 69.920 
30 326 0.152 100 11,551 8.4767 
31 844 0.229 110 3011.0 5.0151 
32 1,274 0.152 100 45,142 33.127 
33 1,115 0.229 90 5768.1 9.6076 
34 615 0.381 110 157.75 0.7273 
35 22 1,408 0.152 100 49,890 36.611 
36 500 0.381 110 128.25 0.59133 
37 300 0.229 90 1552.0 2.5850 

The accuracy of the quadratic approximation was assessed by comparing flows and pressures from a simulation using the approximation to the solution obtained from Hazen–Williams coefficients (Figure 6). Differences between pressures obtained by the models were small, with the largest relative error of 0.62% or 0.24 m at node 14. In the flow solution relative errors had a median of 0% and relative errors in the middle 95% of the data ranged from −12 to 20%. In the remaining 5% of pipes absolute errors ranged from −1.9 to 0.48 L/s. As expected, high relative errors occurred at low flows and had a negligible impact on the overall mass balance.

Figure 6

Comparison of pressures (a) and flows (b) for the 25 node network as using MRE quadratic head loss versus quantities found using the Hazen–Williams coefficients.

Figure 6

Comparison of pressures (a) and flows (b) for the 25 node network as using MRE quadratic head loss versus quantities found using the Hazen–Williams coefficients.

For this small network, the MRE quadratic approximation showed similar computational performance to the conventional Hazen–Williams model. On a Windows laptop with 2.2 GHz processor and 8 Gb RAM, simulations using the quadratic converged in 5 ms (average of 100 runs) and simulations using the Hazen–Williams coefficients converged in 6 ms. Both simulations converged to a tolerance of 10−8 in seven (quadratic) or eight (Hazen–Williams) iterations.

EXNET

To illustrate the applicability of the MRE quadratic approximation on another network, the EXNET water system first mentioned by Farmani et al. (2005) is also considered (Figure 7). This system is based on a real network and is suggested as a benchmark by University of Exeter Centre for Water Systems (2012). Head losses are modeled by the Darcy–Weisbach method. For the purposes of the comparison here, the valves in the network file available online were converted to pipes of length 10 m. A quadratic approximation for the head loss in each pipe was created by choosing at and based on 2 m/s.

Figure 7

The EXNET network modified from Farmani et al. (2005) having 1,894 nodes and 2,469 pipes.

Figure 7

The EXNET network modified from Farmani et al. (2005) having 1,894 nodes and 2,469 pipes.

Flows and pressures calculated using the quadratic approximation compared favorably to the EPANET solution (Figure 8). In the pressure solution, errors ranged from −4.4 to 0.63% and averaged −0.46%. In the flow solution, relative errors had a median of zero and relative errors in the middle 95% of the data ranged from −11 to 3.6%. In the remaining 5% of pipes, flow rates ranged from −15 to 71 L/s and absolute errors ranged from −9.7 to 1.3 L/s. As seen in the figure, these errors did not affect the overall mass balance.

Figure 8

Comparison of pressures (a) and flows (b) for the EXNET system as calculated by using quadratic head loss and Darcy–Weisbach head loss.

Figure 8

Comparison of pressures (a) and flows (b) for the EXNET system as calculated by using quadratic head loss and Darcy–Weisbach head loss.

From a computational perspective, the quadratic approximation showed similar performance to the Darcy–Weisbach model. One hundred simulations converged in an average 15 s for the quadratic model and 17 s for the Darcy–Weisbach model. The quadratic model required 10 iterations and the Darcy–Weisbach model required 11 iterations to converge to a tolerance of 10−8. Profiling of the code indicated that the difference in time was attributable both to the computation of the head loss terms and to the additional iteration for the Darcy–Weisbach model. These results suggest benefits of the quadratic approximation reside not in the computational speed but in the functional form.

CONCLUSIONS

This paper has proposed a quadratic approximation for pipe head loss curves which minimizes the relative approximation error. The approximation applies to pipes and networks parameterized by either the Darcy–Weisbach or Hazen–Williams models. By focusing on the relative rather than absolute error, the approximation stays consistent over the range of possible flows. The accuracy of the approximation depends on the range of pipe flow rates and the pipe roughness. Because rougher pipes develop fully turbulent flow at lower Reynolds numbers, the friction factor varies over a smaller range and a quadratic approximation is more accurate. This result means that the approximation is especially well suited to older water systems where pipes tend to be rougher.

Evaluations of the MRE quadratic approximation proposed here included comparison with other approximations and use in a simulation model. For the Darcy–Weisbach parameterization, MRE quadratic approximations were compared to other essentially quadratic forms. For rough pipes, all approximations show good agreement with the theoretical curve. As pipes become smoother, the MRE quadratic shows closer agreement with the ideal curve. The MRE quadratic approximation was also implemented in a network simulation model and tested on two networks from the literature. Results showed that the MRE quadratic yielded a flow and pressure distribution consistent with results obtained by simulating with the original equations.

In the networks simulated here, the MRE quadratic provided a slight but insignificant improvement in simulation speed compared to the classical equations. Speedup due to the MRE quadratic may become more pronounced for larger networks, especially where thousands of model runs are used to explore the search space of a combinatorial problem such as network expansion or pump scheduling.

Pipe friction models based on quadratic forms open up several new avenues for the theory and practice of water network analysis. Some operational and planning problems can now be formulated as polynomial optimization problems. When additionally exploiting the sparse network structure of an urban water network, this may allow advantage to be taken of sparse semi-definite optimization techniques to approximate the globally optimal solutions of these operational and planning problems.

REFERENCES

REFERENCES
Alperovits
E.
Shamir
U.
1977
Design of optimal water distribution systems
.
Water Resour. Res.
13
(
6
),
885
900
.
American Water Works Association
2004
Steel Pipe: A Guide for Design and Installation
,
4th edn
.
AWWA Manual M11
,
Denver, CO, USA
.
Bragalli
C.
Ambrosio
C. D.
Lee
J.
Lodi
A.
Toth
P.
2011
On the optimal design of water distribution networks: a practical MINLP approach
.
Optim. Eng.
13
(
2
),
219
246
.
Broad
D.
Maier
H.
Dandy
G.
2010
Optimal operation of complex water distribution systems using metamodels
.
J. Water Resour. Plann. Manage.
136
(
4
),
433
443
.
Burgschweiger
J.
Gndig
B.
Steinbach
M.
2009
Optimization models for operative planning in drinking water networks
.
Optim. Eng.
10
,
43
73
.
Christensen
B. A.
Locher
F. A.
Swamee
P. K.
2000
Limitations and proper use of the Hazen–Williams equation
.
J. Hydraul. Eng.
126
(
2
),
167
170
.
Farmani
R.
Savic
D. A.
Walters
G. A.
2005
Evolutionary multi-objective optimization in water distribution network design
.
Eng. Optim.
37
(
2
),
167
183
.
Franzini
J. B.
Finnemore
E. J.
1997
Fluid Mechanics with Engineering Applications
,
9th edn
.
McGraw-Hill
,
New York, USA
.
Gleixner
A.
Held
H.
Huang
W.
Vigerske
S.
2012
Towards globally optimal operation of water supply networks
.
Numer. Algebra Control Optim.
2
(
4
),
695
711
.
Hutton
C.
Kapelan
Z.
Vamvakeridou-Lyroudia
L.
Savic
D.
2014
Dealing with uncertainty in water distribution system models: a framework for real-time modeling and data assimilation
.
J. Water Resour. Plann. Manage.
140
,
169
183
.
Jowitt
P. W.
Xu
C.
1990
Optimal valve control in water-distribution networks
.
J. Water Resour. Plann. Manage.
116
(
4
),
455
472
.
Moody
L. F.
1944
Friction factors for pipe flow
.
Trans. ASME
66
,
671
678
.
Nicolini
M.
Zovatto
L.
2009
Optimal location and control of pressure reducing valves in water networks
.
J. Water Resour. Plann. Manage.
135
(
3
),
178
187
.
Reis
L. F. R.
Porto
R. M.
Chaudhry
F. H.
1997
Optimal location of control valves in pipe networks by genetic algorithm
.
J. Water Resour. Plann. Manage.
123
(
6
),
317
326
.
Rossman
L. A.
2000
EPANET 2 User's Manual
.
US EPA
,
Cincinnati, OH
.
Savic
D. A.
Walters
G. A.
1997
Genetic algorithms for least-cost design of water distribution networks
.
J. Water Resour. Plann. Manage.
123
(
2
),
67
77
.
Sterling
M.
Bargiela
A.
1984
Leakage reduction by optimised control of valves in water networks
.
Trans. Inst. Meas. Control
6
(
6
),
293
298
.
Sterling
M.
Coulbeck
B.
1975
A dynamic programming solution to optimization of pumping costs
.
Proc. Inst. Civil Eng.
59
(
2
),
813
818
.
Swamee
P. K.
Jain
A. K.
1976
Explicit equations for pipe-flow problems
.
J. Hydraul. Div. Am. Soc. Civ. Eng.
102
(
5
),
657
664
.
Todini
E.
Pilati
S.
1988
A gradient method for the analysis of pipe networks
. In:
Computer Applications in Water Supply – Systems Analysis and Simulation
, Vol. 1 (
Coulbeck
B.
Orr
C.-H.
, eds).
Research Studies Press Ltd
,
Letchworth, Hertfordshire, UK
.
Vairavamoorthy
K.
Lumbers
J.
1998
Leakage reduction in water distribution systems: optimal valve control
.
J. Hydraul. Eng.
124
(
11
),
1146
1154
.
White
F. M.
1999
Fluid Mechanics
,
4th edn
.
McGraw-Hill
,
New York, USA
.
Yoo
D. H.
Singh
V. P.
2005
Two methods for the computation of commercial pipe friction factors
.
J. Hydraul. Eng.
131
(
8
),
694
704
.