1. Introduction
The Colebrook function is to date the most used relation in engineering practice for evaluation of flow friction in pipes [
1]. It is given in implicit form
, Equation (1):
In Equation (1), is introduced because of linearization of the unknown flow friction factor λ, Re is the Reynolds number, and is the relative roughness of inner pipe surface (all quantities are dimensionless). Because the Colebrook function in the examined iterative methods sometimes needs to be evaluated in two or three points, y and z are used in the same meaning as x, where they are dimensionless parameters that depend on the friction factor. Practical domain of applicability is for the Reynolds number, and for the relative roughness of the inner pipe surface, .
The Colebrook equation is empirical [
2] and hence its accuracy can be disputed (experiment of Nikuradse or recent experiments from Oregon or Princeton research groups [
3]), but anyway, it is widely accepted in engineering practice. It is based on an experiment conducted by Colebrook and White with the flow of air through a set of pipes with different inner roughness [
2]. The turbulent part of the Moody diagram [
4,
5] is based on the Colebrook equation. The exact solution to this equation does not exist with the exception of those through the Lambert
W-function [
6,
7,
8], but anyway, the Lambert
W-function can be evaluated only approximately [
9,
10,
11]. Many different explicit approximations of the Colebrook equation exist in the form
, and they have different degrees of accuracy and complexity [
12,
13,
14]. Many of the approximations are based on internal iterative cycles [
15,
16,
17], and therefore, it is better to use more accurate iterative procedures if they require only few iterative steps. Calculation of flow through complex networks of pipes, such as for water or gas distribution, requires multiple evaluations of flow friction factor [
18,
19,
20,
21,
22,
23]. In general, the less number of iterations is required, the solution is more efficient with a decreased burden for computers [
24,
25,
26,
27,
28]. Moreover, the recent approach, based on Padé polynomials, shows how the computational burden can be minimized with the same number of used iterations [
29].
The most used iterative methods for solving the Colebrook equation are the Newton–Raphson method or its simplified version, the fixed point methods [
30,
31]. On the other hand, various iterative procedures for solving a single nonlinear equation are available and here are tested, in total, 23 different methods: (i) Fixed-point [
30,
31]; (ii) Newton–Raphson [
32,
33]; Hansen–Patrick group of method described previously [
34] such as: (iii) Halley [
35]; (iv) Euler–Chebyshev [
36,
37]; and (v) Basto–Semiao–Calheiros method [
38]; (vi) Super Halley [
39]; (vii) Ostrowski method of King family [
36,
37,
40]; (viii) Kung–Traub [
41]; (ix) Maheshwari [
42]; (x) Hermite interpolation [
43] based on Jarratt method [
44]; (xi) Khattri and Babajee [
45]; (xii) Murakami [
46]; (xiii) Neta [
47]; (xiv) Chun-Neta [
48] based on Kung–Traub method [
41]; (xv) Wang–Liu [
49] based on Hermit interpolation [
43]; (xvi) Bi–Ren–Wu [
50]; (xvii) Jain [
51] based on Steffensen method [
52,
53,
54,
55]; (xviii) Sharma–Arora [
56]; (xix) Džunić-Petković-Petković [
57,
58,
59]; (xx) Neta–Johnson [
60] based on Jarrat method [
44]; (xxi) Cordero et al. [
54,
55]; (xxii) Sharma–Sharma [
61]; and (xxiii) Sharma–Guha–Gupta method [
62]. All methods [
63] are thoroughly tested within the domain of applicability of the Colebrook equation.
Details about preparation of data for analysis of the iterative methods are given in
Section 2, while the methods are shown in
Section 3 along with numerical examples. Discussion and analysis are given in
Section 4, while concluding remarks and recommendations are in
Section 5.
2. Preparation of Data for Analysis
Calculation of complex networks of pipes for distribution of water or natural gas is a computationally demanding task, where flow friction needs to be evaluated a few million times, which can cause a significant burden for computers. Moreover, a fast but reliable approximation of pipeline hydraulics is needed also for probabilistic risk assessment of pipeline networks. A large number of network simulations of random component failures and their combinations must be automatically evaluated and statistically analyzed in this case [
20,
22]. These tasks can cause significant a burden for computing, which even powerful computer resources cannot easily deal with. In general logarithmic functions noninteger power terms are extremely demanding, and therefore their use should be discouraged [
24,
25,
26,
27]. Both are used frequently in explicit approximations, but noninteger power terms occur very rarely in iterative procedures. So, the main goal of this paper is to find an optimal iterative solution for the Colebrook equation, which leads to the accurate solution within the domain of applicability of the Colebrook equation using the least possible number of iterations. The Colebrook equation is based on logarithmic law and therefore it is hardly possible that the logarithmic function can be eliminated entirely from the iterative procedure [
29].
To maintain compatibility in this analysis with the study of Brkić [
13], which evaluates accuracy of explicit approximations of the Colebrook equation, all iterative methods are tested in 740 uniformly distributed points over the domain of applicability of the Colebrook equation. Examples from five points of 740 in total are shown (five triplets
), of which three are randomly selected, and two are from the most problematic zones which need a slightly higher number of iterations to reach the demanded level of accuracy: (1)
,
; (2)
,
; (3)
,
; (4)
,
; and (5)
,
. The corresponding final solutions are: (1)
; (2)
; (3)
; (4)
; and (5)
. Examples 2 and 5 are those from the more problematic zones of the domain, which in general require more iterative steps to reach the desired accuracy.
The initial starting point for the examined iterative procedures is set to;
. This value is chosen after numerous tests over the domain and it is from the problematic zone, which needs, in general, more iterative cycles to reach the final accurate solution. The fixed initial starting point for the testing is a more suitable option, because the iterative procedures in that way can be compared in a better way. Although the starting point
for the iterative procedure can be chosen using different formulas [
29], numerous tests showed that any fixed value within the domain of applicability of the Colebrook equation leads to the final accurate solution without significant variation of the required number of iterations. From the examined literature, cases when the iterative procedure diverges, oscillates, or converges outside the domain of applicability of the Colebrook equation are not reported.
For the purpose of all examined iterative methods, the Colebrook equation should be in the appropriate form as in Equation (2); for point
;
, where
is the functional symbol for the Colebrook equation. Sometimes it need to be evaluated in additional points
and
. The first
and the second derivatives
are also needed for some methods, but in most cases only in point
. Additional symbols specific for the certain method are explained in
Section 3.
In Equation (2), .
In
Section 3, the presented iterative methods are listed in general from the simplest to the most complex [
63]; (
Section 3.1) One log-call per iteration: (
Section 3.1.1); Fixed-point [
30,
31]; (
Section 3.1.2) Newton–Raphson [
32,
33]; (
Section 3.1.3) Hansen–Patrick [
34]: (
Section 3.1.3a) Halley [
35], (
Section 3.1.3b) Euler–Chebyshev [
36,
37], (
Section 3.1.3c) Basto–Semiao–Calheiros method [
38]; (
Section 3.1.4) Super Halley [
39]; (
Section 3.1.5) Murakami [
46]; (
Section 3.2) Two log-calls per iteration (two-point methods): (
Section 3.2.1) Ostrowski method of the King family [
36,
37,
40]; (
Section 3.2.2) Kung–Traub [
41]; (
Section 3.2.3) Maheshwari [
42]; (
Section 3.2.4) Khattri and Babajee [
45]; (
Section 3.2.5) Hermite interpolation [
43] based on Jarratt method [
44]; (
Section 3.2.6) Wang–Liu method [
49] based on Hermit interpolation [
43]; and finally (
Section 3.3) Three-log-calls per iteration (three-point methods): (
Section 3.3.1) Neta [
47]; (
Section 3.3.2) Chun–Neta [
48] based on Kung–Traub [
41]; (
Section 3.3.3) Džunić–Petković–Petković method [
57,
58,
59]; (
Section 3.3.4) Neta–Johnson [
60] based on Jarrat method [
44]; (
Section 3.3.5) Jain [
51] based on Steffensen method [
52,
53,
54,
55]; (
Section 3.3.6) Bi–Ren–Wu [
50]; (
Section 3.3.7) Cordero et al. [
54,
55]; (
Section 3.3.8) Sharma–Arora [
56]; (
Section 3.3.9) Sharma–Sharma [
61]; and the (
Section 3.3.10) Sharma–Guha–Gupta method [
62].
The complexity of the methods cannot be evaluated only based on the number of log-calls per iteration, because many of them require evaluation of derivatives in one or more points. For example, the first derivative in point
is required for (2) Newton–Raphson, but not for (1) fixed-point method. On the other hand, both methods require the same number of iterations to reach the same level of accuracy. Also, in certain cases, computationally costly log-calls in all subsequent iterations can be substituted with inexpensive Padé polynomials [
64], as described by Praks and Brkić [
29], but also in that way, the number of required iterations to reach the final desired accuracy remained unchanged. Therefore, for evaluation of the efficiency of each of the presented iterative methods, number of iterations are counted in each of the tested 740 points within the domain of applicability. Then, the highest value is chosen as the worst possible representative. We have presented five carefully selected numerical examples that were chosen from these 740, in order to illustrate the calculations.
3. Iterative Methods and Numerical Examples
Iterative methods with one log-call per iteration are presented in
Section 3.1, with two log-calls in
Section 3.2, and three-log calls in
Section 3.3. With possible exceptions, they are sorted with increased complexity. As already explained, numerical examples are given for all presented methods.
In the following formulas indexes i and i + 1 refer to the two subsequent iterations.
In the following tests, the iterative procedures stop when the accuracy in respect to eight decimal places is reached or when sign #div0! appears in the meaning division with zero (desired accuracy reached). Final accurate solutions are marked with shading pattern. Moreover, the final accurate solutions that require the highest number of iterations for the observed methods (the worst cases), are highlighted in the framed cells.
3.1. Iterative Methods with One Log-Call Per Iteration
3.1.1. Fixed-Point Method
Fixed-Point method [
30,
31]; Equation (3):
Example 1: | , | , | , | ; |
Example 2: | , | , | , | , |
| | , | ; | |
Example 3: | , | , | , | ; |
Example 4: | | , | , | ; |
Example 5: | , | , | , | , |
| , | , | , | . |
3.1.2. Newton–Raphson Method
Newton–Raphson method [
32,
33]; Equation (4):
gives the fixed-point method [
30,
31].
Example 1: | , | | | ; |
Example 2: | | | , | , |
| , | , | ; | |
Example 3: | , | , | , | ; |
Example 4: | | , | , | ; |
Example 5: | , | , | , | , |
| , | , | , | . |
3.1.3. Hansen–Patrick Method; (a) Halley; (b) Euler–Chebyshev and (c) Basto, Semiao, and Calheiros Methods
Hansen–Patrick Method [
34] is represented here through (a) Halley [
35]; (b) Euler–Chebyshev [
36,
37] and (c) Basto, Semiao, and Calheiros methods [
38]
(a) Halley method [
35]; Equation (5):
Example 1: | | , | , | ; |
Example 2: | , | , | , | , |
| , | , | ; | |
Example 3: | , | , | , | ; |
Example 4: | , | , | , | ; |
Example 5: | , | , | , | , |
| , | , | , | . |
(b) Euler–Chebyshev method [
36,
37]; Equation (6):
Example 1: | , | , | , | ; |
Example 2: | , | , | , | , |
| , | , | ; | |
Example 3: | , | , | , | ; |
Example 4: | , | , | , | ; |
Example 5: | , | , | , | , |
| , | , | , | . |
(c) Basto-Semiao-Calheiros method [
38]; Equation (7):
Example 1: | , | , | , | ; |
Example 2: | , | , | , | , |
| , | , | ; | |
Example 3: | , | , | , | ; |
Example 4: | , | , | , | ; |
Example 5: | , | , | , | , |
| , | , | | . |
3.1.4. Super Halley Method
Super Halley method [
39]; Equation (8):
Auxiliary parameter is introduced as described.
Example 1: | , | , | , | ; |
Example 2: | , | , | , | , |
| , | , | ; | |
Example 3: | , | , | , | , |
| ; | | | |
Example 4: | , | , | , | ; |
Example 5: | , | , | , | , |
| , | , | , | . |
3.1.5. Murakami Method
Murakami [
46] method; Equation (9):
In addition to the point x, the Murakami method requires an additional evaluation of the function in points and , but only for the first derivative, which does not contain the logarithmic function.
Example 1: | , | , | , | , |
| , | , | , | , |
| , | , | , | ; |
Example 2: | , | , | , | , |
| , | , | , | , |
| , | , | ; | |
Example 3: | , | , | , | , |
| , | , | , | , |
| , | , | , | , |
| ; | | | |
Example 4: | , | , | , | , |
| , | , | , | , |
| , | ; | | |
Example 5: | , | , | , | , |
| , | , | , | , |
| , | , | , | . |
3.2. Iterative Methods with Two Log-Calls Per Iteration
3.2.1. Ostrowski Method of the King Family
Ostrowski method of the King Family [
36,
37,
40]; Equation (10):
Example 1: | , | , | ; | |
Example 2: | , | , | , | ; |
Example 3: | , | , | ; | |
Example 4: | , | , | ; | |
Example 5: | , | , | , | , |
| . | | | |
3.2.2. Kung–Traub Method
Kung–Traub method [
41]; Equation (11):
Example 1: | , | , | ; | |
Example 2: | , | , | , | ; |
Example 3: | , | , | | |
Example 4: | , | , | | |
Example 5: | , | , | , | |
| . | | | |
3.2.3. Maheshwari Method
Maheshwari method [
42]; Equation (12):
Example 1: | , | , | ; | |
Example 2: | , | , | , | ; |
Example 3: | , | , | ; | |
Example 4: | , | , | ; | |
Example 5: | , | , | , | , |
| . | | | |
3.2.4. Khattri and Babajee Method
Khattri and Babajee method [
45]; Equation (13):
Example 1: | , | , | ; | |
Example 2: | , | , | , | ; |
Example 3: | , | , | ; | |
Example 4: | , | , | ; | |
Example 5: | , | , | , | , |
| . | | | |
3.2.5. Hermite Interpolation Based on Jarratt Method
In Equation (14), is Hermite interpolation polynomial. This method requires an evaluation of the first derivatives in points , , and , but the function needs to be evaluated only in points and .
Hermite Interpolation [
43] based on Jarratt Method [
44]; Equation (14):
Example 1: | , | , | , | ; |
Example 2: | , | , | , | , |
| ; | | | |
Example 3: | , | , | , | ; |
Example 4: | , | , | ; | |
Example 5: | , | , | , | , |
| . | | | |
3.2.6. Wang–Liu Method Based on Hermit Interpolation
In Equation (15), is Hermite interpolation polynomial, which is used in the same form as in Equation (14). This method also requires an evaluation of the first derivatives in points , , and , but the function needs to be evaluated only in points and .
Wang–Liu method [
49] based on Hermit Interpolation [
43]; Equation (15):
Example 1: | , | , | , | ; |
Example 2: | , | , | , | |
| , | , | ; | |
Example 3: | , | , | , | ; |
Example 4: | , | , | , | ; |
Example 5: | , | , | , | , |
| , | , | , | . |
3.3. Iterative Methods with Three Log-Calls Per Iteration
3.3.1. Neta Method
Neta method [
47]; Equation (16):
Example 1: | , | | |
Example 2: | , | , | ; |
Example 3: | , | ; | |
Example 4: | , | ; | |
Example 5: | , | , | . |
3.3.2. Chun–Neta Method
Chun–Neta method [
48] based on Kung–Traub [
41]; Equation (17):
Example 1: | , | | | |
Example 2: | , | , | ; | |
Example 3: | , | | | |
Example 4: | , | | | |
Example 5: | , | , | , | . |
3.3.3. Džunić–Petković–Petković Method
Džunić–Petković–Petković method [
57,
58,
59]; Equation (18):
Example 1: | , | ; | |
Example 2: | , | , | ; |
Example 3: | , | ; | |
Example 4: | , | ; | |
Example 5: | , | , | . |
3.3.4. Neta–Johnson Method Based on Jarrat Method
Neta–Johnson method [
60] based on Jarrat method [
44]; Equation (19):
Auxiliary parameter is introduced as described.
Example 1: | , | , | , | , |
| ; | | | |
Example 2: | , | , | , | |
| 4.928634344, | | , | , |
| ; | | | |
Example 3: | , | | , | |
Example 4: | , | , | , | ; |
Example 5: | , | , | , | , |
| , | 4.222050354, | , | , |
| 4.222041013, | , | 4.222041030, | . |
3.3.5. Jain Method Based on Steffensen Method
Jain method [
51] based on Steffensen Method [
52,
53,
54,
55]; Equation (20):
Jain method is a three-point method with evaluation of the function
in points
,
, but also in point
. Serghides explicit approximation of the Colebrook equation is based on Steffensen method [
16,
65].
Example 1: | , | 5.274511499; | |
Example 2: | , | , | ; |
Example 3: | , | ; | |
Example 4: | , | ; | |
Example 5: | , | , | . |
3.3.6. Bi–Ren–Wu Method
Bi–Ren–Wu method [
50]; Equation (21):
Two-parameter function ◊ is introduced as defined.
Example 1: | , | , | ; | |
Example 2: | , | , | ; | |
Example 3: | , | ; | | |
Example 4: | , | , | ; | |
Example 5: | , | , | , | . |
3.3.7. Cordero et al. Method
Cordero et al. method [
54,
55]; Equation (22):
Two-parameter and three-parameter functions ◊ are introduced as defined.
Example 1: | , | , | ; | |
Example 2: | , | , | , | ; |
Example 3: | , | , | ; | |
Example 4: | , | , | ; | |
Example 5: | , | , | , | , |
| . | | | |
3.3.8. Sharma–Arora Method
Sharma–Arora method [
56]; Equation (23):
Two-parameter function ◊ is introduced as defined.
Example 1: | , | ; | |
Example 2: | , | , | ; |
Example 3: | , | ; | |
Example 4: | , | ; | |
Example 5: | , | , | . |
3.3.9. Sharma–Sharma Method
Sharma–Sharma method [
61]; Equation (24):
Two-parameter function ◊ is introduced as defined.
Example 1: | , | ; | |
Example 2: | , | , | ; |
Example 3: | , | ; | |
Example 4: | , | ; | |
Example 5: | , | , | . |
3.3.10. Sharma–Guha–Gupta Method
Sharma–Guha–Gupta method [
62]; Equation (25):
Two-parameter function ◊ is introduced as defined, as well as auxiliary parameters , and .
Example 1: | , | ; | |
Example 2: | , | , | ; |
Example 3: | , | ; | |
Example 4: | , | ; | |
Example 5: | , | , | . |
4. Summary—Discussion and Analysis
To summarize, the highest required number of iterations to reach the final accurate solution of the Colebrook equation in respect to eight decimal places, for the examined methods are:
- 3.1
One log-call per iteration (one-point methods)
- 3.1.1
Fixed-point [
30,
31]; Equation (3): 7 iterations
- 3.1.2
Newton–Raphson [
32,
33]; Equation (4): 7 iterations
- 3.1.3
- 3.1.3a
Halley [
35], Equation (5): 7 iterations
- 3.1.3b
Euler–Chebyshev [
36,
37], Equation (6): 7 iterations
- 3.1.3c
Basto–Semiao–Calheiros method [
38]; Equation (7): 7 iterations
- 3.1.4
Super Halley [
39]; Equation (8): 7 iterations
- 3.1.5
Murakami [
46]; Equation (9): 12 iterations
- 3.2
Two log-calls per iteration (two-point methods):
- 3.2.1
Ostrowski method of the King family [
36,
37,
40]; Equation (10): 4 iterations
- 3.2.2
Kung–Traub [
41]; Equation (11): 4 iterations
- 3.2.3
Maheshwari [
42]; Equation (12): 4 iterations
- 3.2.4
Khattri and Babajee [
45]; Equation (13): 4 iterations
- 3.2.5
Hermite interpolation [
43] based on Jarratt method [
44]; Equation (14): 4 iterations
- 3.2.6
Wang–Liu method [
49] based on Hermit interpolation [
43]; Equation (15): 7 iterations
- 3.3
Three log-calls per iteration (three-point methods):
- 3.3.1
Neta [
47]; Equation (16): 2 iterations
- 3.3.2
Chun–Neta [
48] based on Kung–Traub [
41]; Equation (17): 2 iterations (3 in rare cases)
- 3.3.3
Džunić–Petković–Petković method [
57,
58,
59]; Equation (18): 2 iterations
- 3.3.4
Neta–Johnson [
60] based on Jarrat method [
44]; Equation (19): 11 iterations
- 3.3.5
Jain [
51] based on Steffensen method [
52,
53,
54,
55]; Equation (20): 2 iterations
- 3.3.6
Bi–Ren–Wu [
50]; Equation (21): 3 iterations
- 3.3.7
Cordero et al. [
54,
55]; Equation (22): 4 iterations
- 3.3.8
Sharma–Arora [
56]; Equation (23): 2 iterations
- 3.3.9
Sharma–Sharma [
61]; Equation (24): 2 iterations
- 3.3.10
Sharma–Guha–Gupta method [
62]; Equation (25): 2 iterations
After the conducted analysis the following methods should not be used for the Colebrook equation: (
Section 3.1.5) Murakami [
46]; Equation (9), (
Section 3.2.6) Wang–Liu method [
49] based on Hermit interpolation [
43]; Equation (15), Neta–Johnson [
60] based on Jarrat method [
44]; Equation (19), and Cordero et al. [
54,
55]; Equation (22).
Most one-point methods that require one log-call per iteration need in our case seven iterations to reach the final solution. Among them as the simplest, (
Section 3.1.1) Fixed-point [
30,
31]; Equation (3) can be recommended. Following procedure from Praks and Brkić [
29], only one log-call is required in respect to the whole procedure and in particular only in the first iteration, whereas in all subsequent iterations log-calls are replaced by computationally inexpensive Padé polynomials.
Two log-calls iterative methods (two-point methods) require up to four iterations to reach the final solution. Those simpler can be used: (
Section 3.2.1) Ostrowski method of the King family [
36,
37,
40]; Equation (10), (
Section 3.2.2) Kung–Traub [
41]; Equation (11), and (
Section 3.2.3) Maheshwari [
42]; Equation (12).
Three log-calls methods (three-point methods) are the most powerful among the presented procedures. The simplest, but accurate, can be recommended for use: (
Section 3.3.1) Neta [
47]; Equation (16), (
Section 3.3.2) Chun–Neta [
48] based on Kung–Traub [
41]; Equation (17), (
Section 3.3.3) Džunić-Petković-Petković method [
57,
58,
59]; Equation (18), and (
Section 3.3.5) Jain [
51] based on Steffensen method [
52,
53,
54,
55]; Equation (20). Those which are accurate, but have two-parameter functions ◊, can be used. However, for the Colebrook equation they seem to be too complex and their use should be limited to some rare very well elaborated cases and calculations: Sharma–Arora [
56]; Equation (23), (
Section 3.3.9) Sharma–Sharma [
61]; Equation (24), and (
Section 3.3.10) Sharma–Guha–Gupta method [
62]; Equation (25).
Moreover, for three-point methods, choosing a different starting point
compared with that selected for the shown numerical validation does not alter the number of required iterations significantly [
29,
66]. Thus, the most powerful three-point methods remain recommended even if we changed the initial starting point.
It should be noted that the Colebrook equation is empirical, and that today also many alternative equations exist. These equations are based on various experiments, where the most known are from Princeton or Oregon facilities [
3,
67,
68]. The Princeton facility uses compressed air, while Oregon uses helium, oxygen, nitrogen, carbon dioxide, and sulfur hexafluoride. The Princeton device weights approximately 25 tons, while the Oregon device weighs approximately 0.7 kg.
5. Conclusions
An iterative solution of the Colebrook equation for flow friction in its domain of applicability represents a challenging task, in order to balance speed and accuracy of the iterative process. The fixed-point iterative methods are commonly used, but they demand up to seven iterations to reach the final accurate solution [
30,
31]. On the other hand, numerous explicit approximations with different degrees of accuracy are available [
12,
13,
14]. Similar to the, as here presented, two-point and three-point iterative procedures are very accurate approximations of the Colebrook equation that also contain few internal iterative steps [
15,
16,
17]. Herein we tested and proposed three-point iterative methods that are very accurate, and therefore can be used in a hybrid way: (1) as very accurate explicit approximations if the calculation is terminated after the first iteration; or (2) for very fine computing, in which a very high accuracy can be easily reached after one or two more additional iterations. This is of a great significance in computing of flow through miscellaneous pipe-networks of various applications (water, natural gas, air, etc.), in which multiple evaluations of friction factors take place and when certain cases require only rough estimations; while other very accurate computational models need to be compared, repeated, or analyzed in detail [
69,
70,
71].