Iterative Methods for Looped Network Pipeline

Calculation
Dejan Brkić

To cite this version:

Dejan Brkić. Iterative Methods for Looped Network Pipeline Calculation. Water Resources Manage-
ment, Springer Verlag, 2011, 25 (12), pp.2951-2987. <10.1007/s11269-011-9784-3>. <hal-01586551>

HAL Id: hal-01586551

https://hal.archives-ouvertes.fr/hal-01586551
Submitted on 13 Sep 2017

HAL is a multi-disciplinary open access L’archive ouverte pluridisciplinaire HAL, est

archive for the deposit and dissemination of sci- destinée au dépôt et à la diffusion de documents
entific research documents, whether they are pub- scientifiques de niveau recherche, publiés ou non,
lished or not. The documents may come from émanant des établissements d’enseignement et de
teaching and research institutions in France or recherche français ou étrangers, des laboratoires
abroad, or from public or private research centers. publics ou privés.
 Editorial Manager(tm) for Water Resources Management
Manuscript Draft

Manuscript Number: WARM1450R2

Title: Iterative methods for looped network pipeline calculation

Article Type: Special Issue: 5th Dubrovnik Conference

Keywords: Pipeline networks; Waterworks; Natural gas distribution; Calculation methods; Flow;
Hydraulic systems

Corresponding Author: Dr Dejan Brkic, PhD in Petroleum Eng.

Corresponding Author's Institution:

First Author: Dejan Brkic, PhD in Petroleum Eng.

Order of Authors: Dejan Brkic, PhD in Petroleum Eng.

Abstract: Since the value of the hydraulic resistance depends on flow rate, problem of flow distribution
per pipes in a gas or water distributive looped pipelines has to be solved using iterative procedure. A
number of iterative methods for determining of hydraulic solution of pipeline networks, such as, Hardy
Cross, Modified Hardy Cross, Node-Loop method, Modified Node method and M.M. Andrijašev method
are shown in this paper. Convergence properties are compared and discussed using a simple network
with three loops. In a municipal gas pipeline, natural gas can be treated as incompressible fluid. Even
under this circumstance, calculation of water pipelines cannot be literary copied and applied for
calculation of gas pipelines. Some differences in calculations of networks for distribution of these two
fluids, i.e. water apropos natural gas are also noted.
Response to reviewer's comments
Click here to download Response to reviewer's comments: WARM editor.doc

Dear editor,
All your suggestions are now accepted.
Sincerely yours,
Dejan Brkić, PhD
*Manuscript DB
Click here to download Manuscript: Text DB-rev.doc Click here to view linked References

Iterative methods for looped network pipeline calculation

1
2 Dejan Brkić
3
4
5 PhD, Ministry of Science and Technological Development, Strumička 88, 11050 Beograd,
6
7 Serbia
8
9
10 Tel. +381642543668, e-mail: dejanrgf@tesla.rcub.bg.ac.rs
11
12
13
14
15 Abstract:
16
17 Since the value of the hydraulic resistance depends on flow rate, problem of flow distribution
18
19
20 per pipes in a gas or water distributive looped pipelines has to be solved using iterative
21
22 procedure. A number of iterative methods for determining of hydraulic solution of pipeline
23
24
networks, such as, Hardy Cross, Modified Hardy Cross, Node-Loop method, Modified Node
25
26
27 method and M.M. Andrijašev method are shown in this paper. Convergence properties are
28
29 compared and discussed using a simple network with three loops. In a municipal gas pipeline,
30
31
32 natural gas can be treated as incompressible fluid. Even under this circumstance, calculation
33
34 of water pipelines cannot be literary copied and applied for calculation of gas pipelines. Some
35
36
37 differences in calculations of networks for distribution of these two fluids, i.e. water apropos
38
39 natural gas are also noted.
40
41
42
43
44 Keywords: Pipeline networks, Water distribution system, Natural gas distribution system,
45
46 Calculation methods, Flow, Hydraulic pipeline systems
47
48
49
50
51 1. Introduction
52
53
54 A number of iterative methods for determining the hydraulic solution of water and natural gas
55
56 pipeline networks which take ring-like form, such as, Hardy Cross, Modified Hardy Cross
57
58
59 (including Andrijašev method), Node-Loop method and Modified Node method are compared
60
61
62 1
63
64
65
 in this paper. All presented methods assume equilibrium among pressure and friction forces in
1
2 steady and incompressible flow. As a result, they cannot be successfully used in unsteady and
3
4
5 compressible flow calculations with large pressure drops where inertia force is important.
6
7 Problem of flow in pipes and open conduits was one which had been of considerable interest to
8
9
10 engineers for nearly 250 years. Even today, this problem is not solved definitively. The difficulty
11
12 to solve the turbulent flow problems lies in the fact that the friction factor is a complex function
13
14
15
of relative surface roughness and Reynolds number. Precisely, hydraulic resistance depends on
16
17 flow rate and hence flow problem in hydraulic networks has to be solved iteratively. Similar
18
19 situation is with electric current when diode is in circuit. With common resistors in electrical
20
21
22 circuits where the electric resistances do not depend on the value of electric current in the
23
24 conduit, problem is linear and no iterative procedure has to be used. So problem of flow through
25
26
27 a single tube is already complex. Despite of it, very efficient procedures are available for solution
28
29 of flow problem in a complex pipeline such as looped pipeline like waterworks or natural gas
30
31
32
distribution network is. Most of the shown methods are based on solution of the loop equations
33
34 while the node equations are used only as control of accuracy. Node-Loop method is also based
35
36 on solution of the loop equations but the node equations are also used in calculation and not only
37
38
39 for control purposes. Node method and here shown Modified Node method is based on solution
40
41 of the node equations while the loop equations is used for control purposes.
42
43
44
45
46 In this paper, speed of convergence of presented methods are compared and discussed. This is
47
48
49 done for one simple network with three loops, both, for water and gas distribution networks.
50
51 Similar example can be done for the air ventilation systems in the buildings or mines. For
52
53
54
ventilation problem readers can consult paper of Aynsley (1997), Mathews and Köhler (1995),
55
56 Wang and Hartman (1967), etc. While the disposal of water or gas distribution system is in
57
58
59
60
61
62 2
63
64
65
 single plane with certain elevation of pipes, ventilation network is almost always spatial (Brkić
1
2 2009).
3
4
5
6
7 2. Hydraulics frictions and flow rates in pipes
8
9
10 Each pipe is connected to two nodes at its ends. In a pipe network system, pipes are the channels
11
12 used to convey fluid from one location to another. Pipes are sometimes referred to as tubes or
13
14
15
conduits, and as a part of network to as lines, edges, arcs or branches (in the graph theory, special
16
17 kinds of branches are referred to as trees). The physical characteristics of a pipe include the
18
19 length, inside diameter, roughness, and minor loss coefficient. When fluid is conveyed through
20
21
22 the pipe, hydraulic energy is lost due to the friction between the moving fluid and the stationary
23
24 pipe surface. This friction loss is a major energy loss in pipe flow (Farshad et al. 2001, Kumar
25
26
27 2010).
28
29
30
31
32
Various equations were proposed to determinate the head losses due to friction, including the
33
34 Darcy-Weisbach, Fanning, Chezy, Manning, Hazen-Williams and Scobey formulas. These
35
36 equations relate the friction losses to physical characteristics of the pipe and various flow
37
38
39 parameters. The Darcy-Weisbach formula for calculating of friction loss is more accurate than
40
41 the Hazen-Williams. Beside this, the Hazen-Williams relation is only for water flow.
42
43
44
45
46 2.1 Gas flow rates and pressure drops in pipes
47
48
49 When a gas is forced to flow through pipes it expands to a lower pressure and changes its
50
51 density. Flow-rate, i.e. pressure drop equations for condition in gas distribution networks
52
53
54
assumes a constant density of a fluid within the pipes. This assumption applies only to
55
56 incompressible, i.e. for liquids flows such as in water distribution systems for municipalities (or
57
58 any other liquid, like crude oil, etc.). For the small pressure drops in typical gas distribution
59
60
61
62 3
63
64
65
 networks, gas density can be treated as constant, which means that gas can be treated as
1
2 incompressible fluid (Pretorius et al. 2008), but not as liquid flow. Assumption of gas
3
4
5 incompressibility means that it is compressed and forced to convey through pipes, but inside the
6
7 pipeline system pressure drop of already compressed gas is small and hence further changes in
8
9
10 gas density can be neglected. Fact is that gas is actually compressed and hence that volume of
11
12 gas is decreased and then such compressed volume of gas is conveying with constant density
13
14
15
through gas distribution pipeline. So, mass of gas is constant, but volume is decreased while gas
16
17 density is according to this, increased. Operate pressure for distribution gas network is 4·105 Pa
18
19 abs i.e. 3·105 Pa gauge and accordingly volume of gas is decreased four times compared to
20
21
22 volume of gas at normal conditions.
23
24
25
26
27 Inner surface of polyethylene pipes which are almost always used in gas distribution networks
28
29 are practically smooth and hence flow regime in the typical network is hydraulically „smooth‟
30
31
32
(Sukharev et al 2005). For this regime is suitable Renouard‟s equation adjusted for natural gas
33
34 flow (1):
35
36
4810  Qn1.82  L   r
37 Fg  C  p12  p22 
38 Din4.82
39
40
41 (1)
42
43 Regarding to Renouard‟s formula has to be careful since it does not relate pressure drop but
44
45
46 actually difference of the quadratic pressure at the input and the output of pipe. This means that
47
48 C is not actually pressure drop in spite of the same unit of measurement, i.e. same unit is used
49
50
51 for the pressure (Pa). For gas pipeline calculation is very useful fact that when C →0 this
52
53
54 consecutive means that also C→0. Parameter C can be noted as pseudo-pressure drop.
55
56
57
58
59 2.2. Liquid flow in pipes
60
61
62 4
63
64
65
 In Renoard‟s equation adjusted for gas pipelines (1) friction factor is rearranged in the way to
1
2 be expressed using other flow parameters and also using some thermodynamic properties of
3
4
5 natural gas. Using formulation for Darcy friction factor in hydraulically smooth region
6
7 Renouard suggests his equation for liquid flow (2):
8
9
0.172
10  (2)
11 Re 0.18
12
13
14 Then pressure drop (Ekinci and Konak 2009) can be found very easily using Darcy-Weisbach
15
16 equation (3):
17
18
L 8  Q2
19 Fw  p  p1  p2     
20 D5in  2
21
22
23 (3)
24
25 Renouard‟s equation (2) is based on power law. Liquid flow is better fit using some kind of
26
27
28 logarithmic formula like Colebrook‟s (4):
29
30
1  2.51  
31  2  log   (4)
32   Re  3.71  Din 
33
34
35 Colebrook‟s equation is suitable especially for flow through steel pipes (Colebrook 1939).
36
37 Some researchers adopt a modification of the Colebrook equation, using the 2.825 constant
38
39
40 instead of 2.51 especially for gas flow calculation (Haaland 1983, Coelho and Pinho 2007).
41
42
43
44
45 In the case of waterworks, pressures will also be expressed in Pa, not in meter (m)
46
47 equivalents. In that way, clear comparison with gas distribution network can be done.
48
49 Pressures expressed in Pa can be very easily recalculated in heads expressed in meters (m) for
50
51
52 water networks knowing water density. Example network is located in a single-flat area with
53
54 no variation in elevation (otherwise the correction term must be used which is significant for
55
56
57 water network).
58
59
60
61
62 5
63
64
65
 3. Looped pipeline networks for distribution of fluids
1
2 A pipeline network is a collection of elements such as pipes, compressors, pumps, valves,
3
4
5 regulators, tanks, and reservoirs interconnected in a specific way. In this article focus is on
6
7 pipes. The behavior of the network is governed by the specific characteristics of the elements
8
9
10 and how the elements are connected together. Our assumption is that pipes are connected in a
11
12 smooth way, i.e. so called minor hydraulic loses are neglected. Including other elements
13
14
15
different than pipes is a subject of sufficient diversity and complexity to merit a separate
16
17 review.
18
19
20
21
22 The analysis of looped pipeline systems by formal algebraic procedures is very difficult if the
23
24 systems are complex. Electrical models had been used in study of this problem in the time
25
26
27 before advanced computer became available as background to support demandable numerical
28
29 procedures (Mah and Shacham 1978). Here presented methods use successive corrections or
30
31
32
better to say these methods require iterative procedure. The convergence properties of
33
34 presented methods are the main subject of this article. All of such methods can be divided into
35
36 two groups (1) Methods based on solution of the loop equations, and (2) Methods based on
37
38
39 solution of the node equations.
40
41
42
43
44 Most of the methods used commonly in engineering practice belong to the group based on
45
46 solution of the loop equations. Such method presented in this paper is Hardy Cross method
47
48
49 (Cross 1936). Method of Andrijašev is variation of Hardy Cross method (Andrijašev 1964)1.
50
51
1
52 There are some difficulties with citation of Russian methods. Both methods, Lobačev and
53
54
55 Andrijašev, are actually from 1930‟s but original papers from Soviet time are problem to be
56
57 found. Many books in Russian language explain methods of Lobačev and Andrijašev but in
58
59
60 these books reference list is not available.
61
62 6
63
64
65
 Contemporary with Hardy Cross, soviet author V.G. Lobačev (Latišenkov and Lobačev 1956)
1
2 developed very similar method compared to original Hardy Cross method. Modified Hardy
3
4
5 Cross method proposed Epp and Fowler (1970) which considers entire system simultaneously
6
7 is also sort of loop method. Node-Loop method proposed by Wood and Charles (1972) and
8
9
10 later improved by Wood and Rayes (1981) is combination of the loop and node oriented
11
12 methods, but despite of its name is essentially belong to the group of loop methods. Node-
13
14
15
Loop method is based on solution of the loop equations while the node equations are also
16
17 involved in calculation. Only Node method proposed by Shamir and Howard (1968) is real
18
19 representative of node oriented methods. Node method uses idea of Hardy Cross but to solve
20
21
22 node equations instead of loops ones. Improved Node method uses similar idea as Epp and
23
24 Fowler (1970) suggested for the improvement of original Hardy Cross method, but of course
25
26
27 here applied to the Node method.
28
29
30
31
32
Example network with three loops is shown in the figure 1. Before calculation, maximal
33
34 consumption for each node including one of more inlet nodes has to be determined. Pipe
35
36 diameters and node inputs and outputs cannot be changed during the iterative procedure. Goal
37
38
39 is to find final flow distribution for this pipeline system (simulation problem).
40
41
42
43
44 Figure 1. Example of pipeline network
45
46
47
48
49 Pipe lengths and pipes diameters are listed in figure 1. Final flows do not depend on first
50
51 assumed fluid flows per pipes in the case of loop oriented methods or on first assumed
52
53
54
pressure drops or pseudo-pressure drops in the case of node oriented methods. To have a
55
56 better appreciation of the utility of these representations, first will be considered the laws that
57
58 govern flow rates and pressure drops in a pipeline network. These are the counterparts to
59
60
61
62 7
63
64
65
 Kirchhoff‟s laws for electrical circuits and accordingly for hydraulic networks, namely, (i) the
1
2 algebraic sum of flows at each node must be zero; and (ii) the algebraic sum of pressure drops
3
4
5 for water network (i.e. pseudo-pressure drop for gas network) around any cyclic path (loop)
6
7 must be approximately zero. For a hydraulic network, fact that both laws are satisfied
8
9
10 simultaneously means that calculation for steady state for simulation problem is finished.
11
12
13
14
15
3.1 Methods based on solution of the loop equations
16
17 Loops are sometimes referred also to as contours or paths. Note that contours and loops are
18
19 not synonyms in the M.M. Andrijašev method.
20
21
22
23
24 For the loop oriented methods first Kirchhoff‟s law must be satisfied for all nodes in all
25
26
27 iterations. Second Kirchhoff‟s law for each loop must be satisfied with acceptable tolerance at
28
29 the end of the calculation. Or, in other words, final flows are these ones which values are not
30
31
32
changed between two successive iterations (must be satisfied for flow in each pipe).
33
34
35
36 In a loop oriented methods, first initial flow pattern must be chosen to satisfied first
37
38
39 Kirchhoff‟s law. Endless number of flow combinations can satisfy this condition. Someone
40
41 can conclude that it is crucial to start with good initial guesses. But, how does one obtain
42
43
44 good initial guesses? And how sensitive are the methods to the initial guesses? Answer is that
45
46 initial flow pattern does not have significant influence on convergence properties of observed
47
48
49 method. These indicate that the choice of initial flows is not critical and need to be chosen
50
51 only to satisfied first Kirchhoff‟s law for all nodes (Gay and Middleton 1971). It would be
52
53
54
equally satisfactory to generate the initial distribution within the computer program. Initial
55
56 flow pattern for our example network is shown in the figure 2. Loops are also defined in the
57
58 figure 2.
59
60
61
62 8
63
64
65
 1
2 Figure 2. Example network initial parameters prepared for loop oriented calculation
3
4
5
6
7 Results of calculation using original Hardy Cross, Modified Hardy Cross, M.M. Andrijašev
8
9
10 method are not actually flows Q but correction of flows ΔQ calculated for each loop. These
11
12 corrections have to be added algebraically to flows from previous iteration for each pipe
13
14
15
according to specific rules. A pipe common to two loops receives two corrections. The upper
16
17 plus or minus sign (shown in tables 1-8) indicates direction of flow in that conduit in these
18
19 two loops and is obtained from Q for previous iteration. The upper sign is the same as the sign
20
21
22 in front of Q if the flow direction in each loop coincides with the assumed flow direction in
23
24 the particular loop under consideration, and opposite if it does not. The lower sign is copied
25
26
27 from the primary loop for this correction (sign from the loop where this correction is first,
28
29 sign preceding the first iteration from adjacent contour for the conduit taken into
30
31
32
consideration). The rules for sign of corrections ΔQ2 are: (1). the algebraic operation for
33
34 correction ΔQ1 should be the opposite of its sign; i.e. add when the sign is minus. (2). the
35
36 algebraic operation for corrections ΔQ2 should be the opposite of their lower signs when their
37
38
39 upper signs are the same as the sign in front of Q, and as indicated by their lower signs when
40
41 their upper signs are opposite to the sign in front of Q. For details of sing of corrections
42
43
44 consult paper of Brkić (2009) and article Corfield et al (1974) from Gas Engineers Handbook.
45
46
47
48
49 3.1.1 Hardy Cross method for the pipeline calculation
50
51 Hardy Cross, American engineer, developed method in 1936, later named after him (Cross
52
53
54
1936). Hardy Cross calculation for gas pipeline network is shown in table 1 and for water
55
56 network in table 2. Only two iterations will be shown in details. Flow Q 1 calculated in first
57
58 iteration become initial flow Q for second iteration. The plus or minus preceding the flow Q,
59
60
61
62 9
63
64
65
 indicates the direction of the conduit flow for the particular loop (Figure 2). A plus sign
1
2 denotes clockwise flow in the conduit within the loop; a minus sign, counter-clockwise. A
3
4
5 flow correction ΔQ1 as shown in table 1 and 2 is computed for each loop. This correction
6
7 must be subtracted algebraically from the assumed gas flow according to specific rules
8
9
10 explained previously.
11
12
13
14
15
Correction for each loop ΔQ1 is calculated using first derivative of pressure function C for
16
17 pipes defined by Renouard‟s equation for gas networks (5) and derivative of pressure drop Δp
18
19 for pipes defined by Darcy-Weisbach for water networks (6) where flow Q is treated as
20
21
22 variable:
23
24
 4810  Qn1.82  L   r 
 
 
25
C  Din4.82    R Q
1.82
1.82  4810  Qn0.82  L   r
 
26
27 F'   1.82  R  Qn0.82
28 Q Q Q Din4.82

29
30 (5)
31
32
33  8      L  Q2 
 
34
35 F'
p 
   2  D5  
   R  Q  16      L  Q  2  R  Q
2

36 Q Q Q  2  D5
37
38
39 (6)
40
41 Derivative for loops is calculated using assumed loop shown in the figure 2 with no reference
42
43
44 to direction. For the loop I this derivative for gas network (7) and for water network (8) can
45
46 be written:
47
48
49
50
CI Q    C3 Q3   C4 Q4    C7 Q7  I
51
 Q 0.82  L Q 0.82  L Q 0.82  L 
52 FI'    1.82  4810   r   3 4.82 3  4 4.82 4  7 4.82 7 
53 QI QI  D3 D4 D7 I
54
55
56 (7)
57
58
59
60
61
62 10
63
64
65
 pI Q    p3 Q3   p4 Q4    p7 Q7  I 16    3  L3  Q3 4  L 4  Q4 7  L 7  Q7 
FI'    2     
1 QI QI   D53 D54 D57 I
2
3
4
5 (8)
6
7 In presented example loop I begins and ends in node II via pipes 3, 4, and 7.
8
9
10
11
12 Correction of flow (ΔQ1) for each loop can be calculated (9):
13
14
F
15 Q1  
16 F' (9)
17
18
19
20
21
22 Table 1. Hardy Cross calculation for example gas network
23
24
25
26
Table 2. Hardy Cross calculation for example water network
27
28
29
30
31 In the original Hardy Cross method, each loop correction is determined independently of
32
33
34 other loops. Original Hardy Cross method is a sort of Newton–Raphson method but used to
35
36 solve each single loop equation solely, one by one. Hence, Hardy Cross method is also known
37
38
39 as the single contour adjustment method. In matrix form, original Hardy Cross method for the
40
41 example network for water distribution from figure 1 can be noted as (10):
42
 FI   Q3 , Q4 ,  Q7 
43
44 
 0 0 
45
  QI    QI   FI 
46
 
FII Q1 ,  Q 2 ,  Q4 , Q5       (10)
47  0
 QII 
0  x  QII    FII 
48   Q   F 
49  FIII   Q5 , Q6 , Q7 ,  Q8   III   III 
 0 0 
 QIII 
50
51  
52
53 Using numerical values from table 2 for the first iteration of calculation of water network matrix
54
55
56 equation can be written (11):
57
58
59
60
61
62 11
63
64
65
 129435625 0 0   QI  - 15928498
1  0 53910587 0  x  Q    5040086  (11)
2    II   
3  0 0 193795963 QIII   2070510 
4
5
6 Then set of correction for water network calculation in the first iteration is [ΔQI, ΔQII, ΔQIII]T=[-
7
8
9
12306·10-5, 9349·10-5, 1068·10-5]T which is identical as in table 2.
10
11
12
13 3.1.2 Modified Hardy Cross method for the pipeline calculation
14
15
16 Modified Hardy Cross method (somewhere called improved) is also known as the
17
18 simultaneous loop adjustment method. As seen in figure 1, several loops have common pipes,
19
20
21 so corrections to these loops will cause energy losses around more than one loop. In figure 1,
22
23 pipe 4 belongs to two loops (loop I and II), pipe 7 to loop I and III, and finally pipe 5 to II and
24
25
26
III. Modified Hardy Cross method is a sort of Newton–Raphson method used to solve
27
28 unknown flow corrections taking into consideration whole system simultaneously. Epp and
29
30 Fowler (1970) gave idea for this approach. To increase efficiency of Hardy Cross method
31
32
33 zeros from non-diagonal term in matrix equation (10) will be replaced to include influence of
34
35 pipes mutual with adjacent loop (12).
36
37
38  FI   Q3 , Q4 ,  Q7  F   Q4  F  Q  
  I  I 7 
39
  QI   QII   QIII    Q   F 
FII   Q4  FII  Q1 ,  Q2 ,  Q4 , Q5  FII   Q5 
40
41    I
 
I

    x  QII    FII 
42
   Q I    Q II    Q III   Q   F 
FIII  Q7  FIII   Q5  FIII   Q5 , Q6 , Q7 ,  Q8   III   III 
43 
44    
45
  QI   QII   QIII  
46
47 (12)
48
49 Presented matrix is symmetrical; for example (13):
50
51
FI (Q7 ) FIII (Q7 )
52  (13)
53 QIII  QI 
54
55
56
57
58
59
60
61
62 12
63
64
65
 This is because pipe 7 is mutual for two adjacent loops (loop I and loop III). Non-diagonal
1
2 terms have always opposite sign than diagonal. Spatial networks common for ventilation
3
4
5 systems in buildings or mines are exceptions (Brkić 2009).
6
7
8
9
10 Using numerical values from table 2 for the first iteration of calculation of water network matrix
11
12 equation can be written (14):
13
14
15 12943562525 - 4978573 - 70157780  QI  - 159284985
16  - 4978573 53910587 - 429677  x  QII    50400861 
17 
18  - 70157780 - 429677 193795963 QIII   20705102 
19
20
21 (14)
22
23 First two iterations from our example using the modified Hardy Cross method is shown in
24
25
26
table 3 for gas network and in table 4 for water network.
27
28
29
30 Table 3. Calculation after the modified Hardy Cross method for example gas network
31
32
33
34
35 Table 4. Calculation after the modified Hardy Cross method for example water network
36
37
38
39
40 3.1.2 Modified Andrijašev method for the pipeline calculation
41
42
43 Andrijašev method can be used in the formulation as in the original Hardy Cross method and
44
45 as in Modified Hardy Cross method. Here it will be given in notation as improved method
46
47
because this approach shows better convergence performance (example for two iteration of
48
49
50 gas network calculation is shown in table 5 and for water network in table 6). It can be noted
51
52 that some pipes in table 1-4 received only one correction per iteration (for example pipe 3 in
53
54
55 contour I). This means that pipe 3 belongs only to one loop. In method of M.M. Andrijašev
56
57 contours can be defined to include few loops. Thus, contours can be defined in other way and
58
59
60 then each pipe in the network belongs to two contours (Figure 3). Loop is not synonym with
61
62 13
63
64
65
 contour in M.M. Andrijašev method as in Hardy Cross approach. Andrijašev‟s contour will be
1
2 marked with special sign (O).
3
4
5
6
7 Figure 3. Contours for method of M.M. Andrijašev calculation
8
9
10
11
12 Now contour IO (red circuit in figure 3) starting and ending in node I via pipes 4, 5. 6, 8, 3,
13
14
15
contour IIO (green circuit in figure 3) starting and ending in referent node via pipes 1, 6, 8, 7,
16
17 4, 2, and finally contour IIIO (blue circuit in figure 3) starting and ending in referent node via
18
19 pipes 1, 5, 7, 3, 2.
20
21
22
23
24 Matrix formulation of this method for example gas network can be written as (15):
25
26
27  FIo   Q3 , Q4 ,  Q5 , Q6 ,  Q8  FIo   Q4 , Q6 , Q8  FIo  Q3 ,  Q5  
28

   Q  o
 Q 
o
  Q o
 
  QIo   FIo 
FII   Q4 , Q6 , Q8  FII   Q1 ,  Q2 ,  Q4 , Q6 , Q7 , Q8  FII  Q1 , Q2 ,  Q7 
I II III

29  o o o
  o   o
 x  QII    FII 
30 
   Q  o
 Q 
o
  Q o
  Q o   F o 
FIII  Q3 ,  Q5  FIII  Q1 , Q2 ,  Q7  FIII  Q1 ,  Q2 ,  Q3 , Q5 ,  Q7   III   III 
I II III

o o o
31
32 
  
 QIo 
 QIIo 
 QIIIo  

33
34 (15)
35
36 Here has to be very careful because non-diagonal terms are not always negative as in
37
38
39
modified Hardy Cross method (16):
40
40489565894 25862392143 11719977603  QI  - 3277663792
41 o

42  25862392143 52949765748 - 4166255733 x Q o    1492371894 

 (16)
43    II   
44 11719977603 - 4166255733 36379348868 QIIIo  - 2518358643
45
46
47 For example term in first row and second column from (16) is
48
49
25862392143=4298435730+22897756035-1333799622 according to table 5. Same value has
50
51
52 term in second row, first column, etc. Presented matrix is symmetrical. Similar can be done
53
54 for water network (table 6).
55
56
57
58
59 Table 5. Calculation after Modified M.M. Andrijašev method for example gas network
60
61
62 14
63
64
65
 1
2 Table 6. Calculation after Modified M.M. Andrijašev method for example water network
3
4
5
6
7 M.M. Andrijašev method does not improve convergence properties compared with Modified
8
9
10 Hardy Cross method. It only complicates calculation of non-diagonal terms in the matrix
11
12 which is used for calculation. Furthermore, calculation of single-flat pipeline networks using
13
14
15
M.M. Andrijašev method is equally complex as calculation of spatial, multidimensional
16
17 network, using Modified Hardy Cross method (Brkić 2009).
18
19
20
21
22 3.2 Methods based on solution of the node equations
23
24 Nodes are sometimes also referred to as junctions, points or vertices.
25
26
27
28
29 Shamir and Howard (1968) introduced first node oriented method. For the node oriented
30
31
32
methods second Kirchhoff‟s law must be satisfied for each loop in all iterations. First
33
34 Kirchhoff‟s law for each node must be satisfied with acceptable tolerance at the end of the
35
36 calculation. Here will be presented only one method. Improvement of presented method is
37
38
39 done according to the same idea as used to improve original Hardy Cross method but here
40
41 applied to the matrix solution for the node equations.
42
43
44
45
46 3.2.1 Node method
47
48
49 Pipe equations in previous text were expressed as Δp=F(Q) for waterworks, or
50
51 C= p12  p22 =F(Q) for gas networks. These relations can be rewritten in form as Q=f(Δp) for
52
53
54 waterworks or Q=f(C) for gas networks. After that, Renouard‟s equation (1) can be
55
56 rearranged (17):
57
58
59
60
61
62 15
63
64
65
  
1 1
 p12  p22  Din4.82 1.82  C  Din4.82 1.82
1
f g  Q(C )      
2
 4810  L   r   4810  L   r 
(17)
3
4
5 Similar, can be done for Darcy-Weisbach equation (18):
6
7
8
f w  Q(p) 
 p1  p2  D5in   2 
p  D5in   2
9 8  L   8  L   (18)
10
11
12 In Node method, (C) for each pipe in a gas network and (Δp) for each pipe in a water network
13
14
15
has to be assumed, not flows. These assumed functions of pressure must be chosen to satisfy
16
17 second Kirchhoff‟s law (Figure 4). First Kirchhoff‟s law will be fulfilled with demanded
18
19 tolerance at the end of calculation. For gas network, correction of C, noted as ΔC must be
20
21
22 calculated, and for water network, correction of Δp, here noted as ΔΔp. Same algebraic rules as
23
24 for loop oriented methods are valid. Correction ΔC or ΔΔp will be calculated using first
25
26
27 derivative of Q=f(Δp) for waterworks and Q=f(C) for gas networks, where Δp, i.e. C is treated
28
29 as variable (19):
30
31
f
32 x  
33 f' (19)
34
35
36 In previous equation x is Δp or C. Example network adjusted for node oriented calculation is
37
38
39
shown in figure 4 (red letters for gas network and blue for water network).
40
41
42
43 Figure 4. Example of pipeline gas network with three loops adjusted for node oriented
44
45
46 methods
47
48
49
50
51 Calculation of the network using non-improved Node method is not shown here.
52
53
54
55
56 3.2.2 Modified Node method
57
58 Modified Node method is also referred to as the simultaneous node adjustment method. For
59
60
the improvement of Node method, same idea as for the improvement of the Hardy Cross
61
62 16
63
64
65
 method is used. One node must be chosen to be omitted from the calculation shown in table 7
1
2 and 8, because linear independency among node equations is preserved in that way (Mathews
3
4
5 and Köhler 1995, Mah and Shacham 1978). More details about graph theory of networks will
6
7 be explained in further text (section 3.4). Gas network calculation after Modified Node
8
9
10 method is shown in table 7:
11
12
13
14
15
Table 7. Calculation after Modified Node method for example gas network with three loops
16
17
18
19 Corrections in the table 7 and 8 are calculated using matrix equation (20):
20
21
 f1'   x2 ,  x3 ,  x4   f1' x3   f1' x4    x  
  1   2
22 0 0 Q  Q3  Q 4  QI 

 f 2' x3  f   x3 ,  x7 ,  x8   f x7   f x8   x2 
    Q3  Q 7  Q8  QII 
' ' '
23  2 2 0 2
24 
  f 3' x4   f 3' x7  f 3'   x4 ,  x5 ,  x7   f 3'  x5  0 
 x  x3     Q 4  Q5  Q 7  QIII 
   
f 4'   x1 ,  x5 ,  x6 
25   f 4' x5   f 4 x6    4   1
    IV 
'
0 0 x Q Q 5 Q 6 Q
26  
f 5'   x6 ,  x8   x5    Q 6  Q8  QV
 
27  0  f x8 
5
'
0  f 5' x6  
28 (20)
29
30
31 When non-diagonal terms in the first matrix (21) are equalized with zero, results are equal as for
32
33 Node method in basic form (19). For the first iteration, numerical values for the matrix equation
34
35
36 are extracted from table 6 (21):
37
38  3.69 10-10 - 6.21 10-11 - 6.28 10-11 0 0   x1    0.7665
39      
40  - 6.21 10
-11
1.28 10-10 - 4.19 10-11 0 - 2.46 10-11  x2   - 0.9121 
41 - 6.28 10-11 - 4.18 10-11 3.40 10-10 - 2.36 10-10 0  x  x3    2.9302
42  -11     
43  0 0 - 2.36 10 -10
3.54 10 -10
- 3.30 10  x4   - 2.2805 
 - 2.46 10-11 - 3.30 1011 5.76 10-11   x5   0.4417 
44  0 0
45
46
47 (21)
48
49 For the node oriented calculation applied for water networks, pressure function (C) must be
50
51
52 replaced with pressure drop (Δp). Further analogy is clear (Table 8).
53
54
55
56
57 Table 8. Calculation after Modified Node method for example water network with three loops
58
59
60
61
62 17
63
64
65
 3.4 Node-Loop method
1
2 Conditionally, Node-Loop method can be sorted in the group of methods based on solution of
3
4
5 the loop equations according to the previous discussion, but better solution is to be treated as
6
7 combination of the loop and node oriented methods as its name unambiguously suggest.
8
9
10 Node-Loop method is also known as the flow adjustment method.
11
12
13
14
15
Wood and Charles (1972) developed the flow adjustment method by coupling the loop
16
17 equations with the node equations. Wood and Rayes (1981) improved this method. Rather
18
19 than solve for loop corrections, in this method, conservation of energy around a loop is
20
21
22 written directly in the terms of the pipe flow rates. Final result after this method is not flow
23
24 correction, but even better flow itself.
25
26
27
28
29 Pipeline network with three loops and six nodes is used as example in this article. Graph has
30
31
32
X branches (pipes) and Y nodes where in our example, X=8 and Y=6. Graph with Y nodes (in
33
34 our case 6) has Y-1 independent nodes (in our case 5) and X-Y+1 independent loops (in our
35
36 case 3). Tree is a set of connected branches chosen to connect all nodes, but not to make any
37
38
39 closed path (not to form a loop). Branches, which do not belong to a tree, are links (number of
40
41 links are X-Y+1). Loops in the network are formed using pipes from tree and one more
42
43
44 chosen among the link pipes. Number of the loops is determined by number of links. Network
45
46 from example has six nodes and five independent nodes. One node can be omitted from
47
48
49 calculation while no information on the topology in that way will be lost. Rows in the node
50
51 matrix with all node included are not linearly independent. To obtain linear independence any
52
53
54
row of the node matrix has to be omitted. For example, pipe 6 is between node IV and V, and
55
56 reasonable assumption is that if node IV is output node for flow through pipe 6, then node V
57
58 must be input node for flow through this pipe. In our example, node VI will be omitted. First
59
60
61
62 18
63
64
65
 Kirchhoff's law for the initial flow pattern shown in figure 2 can be written using set of
1
2 equations (22):
3
4
5 nodeI ~ QI  Q2  Q3  Q4  0
6
7 nodeII ~ QII  Q3  Q7  Q8  0
8
9
10
nodeIII ~ QIII  Q4  Q5  Q7  0
11
12
nodeIV ~ QIV  Q1  Q5  Q6  0 (22)
13
14
nodeV ~ QV  Q6  Q8  0
15
16
nodeVI ~ QVI  Q1  Q2  0
17
18 Node VI will be omitted from the node matrix to assure linear independency of the rows as
19
20
21 shown in figure 4 (23):
22

 0 1 1 1 0 0 0 0 
23
24
25  0 0  1 0 0 0  1  1
26
 
N   0 0 0  1  1 0  1 0 
27
28 
29   (23)
30
  1 0 0 0  1  1 0 0 
31
32  0 0 0 0 0  1 0  1
33
34
35 First Kirchhoff's law must be fulfilled in all iterations for all nodes. Second Kirchhoff's law
36
37
38 for the initial water flow pattern shown in the figure 1 can be written using set of equations
39
40 (24):
41
42
43 loop I ~  F3  F4  F7  FI
44
45
46
loop II ~ F1  F2  F4  F5  FII
(24)
47
48 loop III ~  F5  F6  F7  F8  FIII
49
50
51 Second Kirchhoff's law for the initial flow pattern shown in the figure 1 also can be noted in
52
53 matrix form using loop matrix (25):
54
55
56
57
58
59
60
61
62 19
63
64
65
 1  0 0 1 1 0 0 1 0 
2
3 L   1  1 0  1  1 0 0 0  (25)
 0 0 0 0  1  1  1  1
4
5
6
7
8 Second Kirchhoff's law must be fulfilled for all loops at the end of calculation with demanded
9
10 accuracy (i.e. FI→0, FII→0 and FIII→0).
11
12
13
14
15 In Node-Loop method these two matrices become one with some modifications. The nodes
16
17
18 and the loops equations shown in previous text here will be united in one coherent system by
19
20 coupling these two set of equations. To introduce matrix calculation, the node-loop matrix
21
22 [NL], matrix of calculated flow in observed iteration [Q], and [V] matrix will be defined (26):
23
24
25  QI 
26
 
27
  QII 
28
29  QIII 
30  
 QIV 
V   
31
32
QV  (26)
 
 
33
34
  FI   F3  F3'  F4  F4'  F7  F7' 
 
35
36   F  F  F'  F  F'  F  F'  F  F' 
 
 
II 1 1 2 2 4 4 5 5
37
38  FIII   F5  F5  F6  F6  F7  F7  F8  F8'
' ' ' 
39  
40
41 Sign minus preceding some of the flows Q in matrix [V] means that this particular Q is not
42
43
44 consumption (sing minus represent inlet of fluid). Node-Loop matrix [NL] can be defined
45
46 using node matrix, loop matrix and first derivative of Renouard‟s function for gas pipes or of
47
48
49 Darcy-Weisbach for water networks, as follows (27):
50
51
52
53
54
55
56
57
58
59
60
61
62 20
63
64
65
  0 1 1 1 0 0 0 0 
1  0 0 1 0 0 0 1 1 
2  
3  0 0 0 1 1 0 1 0 
4  
1 0 0 0 1 1 0 0
NL   0 
5
6 0 0 0 0 1 0 1 
7  
8  0 0  1  F3'  1 F4' 0 0  1 F7' 0 
9  1  F '  1 F2' 0  1 F4'  1  F5' 0 0 0 
10  1

11  0 0 0 0  1 F5'  1  F6'  1 F7'  1 F8' 
12 
13
14 (27)
15
16
17 Further, vector [Q] of the unknown flows can be calculated in the first iteration (28).
18
19
 Q1 
Q 
20
21
22  2
23 Q3 
24  
Q
25
26
Q   4   invNL xV  (28)
Q5
27  
28 Q6 
29
30
Q 
31  7
32 Q8 
33
34
35
Possible sign minus in a front of flow Q in the matrix [Q] means that calculated flow direction
36
37 is opposite compared to shown one in the previous iteration (or in the figure 1 in our case for
38
39 the first calculated values of flows compared with initials flow pattern). If all values of
40
41
42 pressure drops sums calculated after (24) are not approximate zero with reasonable accuracy,
43
44 calculation has to be repeated using values calculated in previous iteration. At the end of
45
46
47 calculation calculated set of flows [Q] stays unchanged.
48
49
50
51
52
Here will be used values from table 1 for gas network (29) and from table 2 for water network
53
54 (30) as example. These values are valid for the first iteration. First five rows (first matrix) are
55
56 from node equation, and next three is from loop equation but multiplied with first derivate
57
58
59 marked in tables as F‟:
60
61
62 21
63
64
65
  0 1 1 1 0 0 0 0   Q1   0.055 
 0 0  1 0 0 0  1  1  Q   - 0.277 
1    2  
2  0 0 0 1 1 0 1 0  Q3   0.361 
     
3  1 0 0 0  1  1 0 0 Q
 x 4    0.222 
4  0 0 0 0 0 1 0 1  Q5   0.194 
     
5  0 0 - 11839776055 1333799622 0 0 - 14293015047 0  Q6  - 2988241676
8812326713 - 1314432601 0 - 1333799622 119798452.1 0 0 0  Q   923187587 
6    7  
7  0 0 0 0 - 119798452.1 4298435730 14293015047 - 22897756035 Q8   300557365 
8 (29)
9
10  0 1 1 1 0 0 0 0   Q1   0.055 
11  0 0  1 0 0 0  1  1  Q   - 0.277 
12    2  
 0 0 0 1 1 0 1 0  Q3   0.361 
13      
14  1 0 0 0 1 1 0 0  x Q4    0.222 
 0 0 0 0 0 1 0 1  Q5   0.194 
15      
16  0 0 - 54299272 4978573 0 0 - 70157780 0  Q6  - 15928498
17 426229242 - 5879434 0 - 4978573 429677 0 0 0  Q   5040086 
   7  
18  0 0 0 0 - 429677 16717415 70157780 - 106491089  Q8   2070510 
19
20
(30)
21 For gas network last three rows are calculated as follows using values from table 1,
22
23
24 for loop I; -2988241676=
25
26 =-3644197165+(-0.194·11839776055+0.027·1333799622+(-0.305)·14293015047)
27
28
29
for loop II; 923187587=1125838521+
30
31 +(0.277·8812326713+(-0.277·1314432601)+(-0.027)·1333799622)+0.027·119798452)
32
33 for loop III; 300557365=366533372+
34
35
36 +(-0.027·119798452+0.027·4298435730+0.305·14293015047+(-0.166·22897756035))
37
38
39
40
41 After first iteration for gas network vector of flows is [0.198409265, 0.357146291,
42
43 0.043307855, 0.25828288, -0.094469817, 0.07065686, 0.197298048, 0.123787585]T. Minus
44
45
46 in front of flow in pipe 5 means: change assumed flow direction from previous iteration. After
47
48 first iteration for water network vector of flows is [0.197719798, 0.357835758, 0.052496097,
49
50
0.249784106, -0.092806697, 0.068304272, 0.204133702, 0.126140172]T. Flows are
51
52
53 expressed in m3/s.
54
55
56
57
58 Excellent book for waterworks calculation by Boulos et al. (2006) can be recommended for
59
60 further reading. In this book, authors instead to omit one node in the node matrix to preserve
61
62 22
63
64
65
 linear independency of rows in this matrix introduce one pseudo-loop in loop matrix. This
1
2 procedure is not practical because at least two nodes with equal pressure must be found in the
3
4
5 network. This is not always possible. Further in that way the node-loop matrix has two additional
6
7 rows which could be avoided. Mathews and Köhler (1995) in his discussion use simplest way,
8
9
10 i.e. they omit one row.
11
12
13
14
15
4. Comparison of solution techniques for looped piping networks
16
17 Final flows are unique after all presented methods, and will be listed in table 9, both for water
18
19 and for gas network.
20
21
22
23
24 Table 9. Final flows for network presented in this paper
25
26
27
28
29 Each method has advantages and shortcomings. Convergence performances will be compared
30
31
32
for all presented methods in figure 5. Note that Modified Node method cannot be compared
33
34 literary because initial values cannot be equalized. In all other methods initial patterns are
35
36 given in the form of flows and equalized, while in Node method initial pattern is in the form
37
38
39 of pressures (better to say function of pressure).
40
41
42
43
44 Figure 5. Comparisons of convergence for presented methods (gas network)
45
46
47
48
49 Best way to compare water and gas distribution network is to compare velocity of gas and
50
51 water through pipes. For water, this can be done using (31):
52
53
4Q
54 vw  (31)
55 Din2  
56
57
58 Velocity of gaseous fluids depends on the pressure in pipe since they are compressible (32):
59
60
61
62 23
63
64
65
 4  pn  Qn
vg  (32)
1 p  Din2  
2
3
4 Velocities for water and gas for calculated flows through the pipes in our example is listed in
5
6
7
table 10:
8
9
10
11 Table 10: Velocities for water and gas for calculated flows from example network
12
13
14
15
16 5. Conclusion
17
18
19 Comparison between analyzed methods was carried out, taking as a criterion the number of
20
21 iteration to achieve final results. Modified Hardy Cross method, Modified M.M. Andrijahshev
22
23
24 method and Node-Loop method have equal performances according to above adopted criterion.
25
26 For more complex networks, using Node-Loop method, number of required iterations is smaller
27
28
even compared with Modified Hardy Cross method. Among these three methods, Node-Loop
29
30
31 method is superior because it does not require complex numerical scheme for algebraic addition
32
33 of corrections in each of iterations. In Node-Loop method, final results of each of the iterations
34
35
36 are flows directly and not correction of flows. Modified Andrijašev method is more complex
37
38 compared with Modified Hardy Cross method but without any improvement in the properties of
39
40
41 the convergence. Node method has the worst performance of convergence, but this method is
42
43 different in its approach compared with all other methods shown in this paper. Node method
44
45
46
cannot be rejected easily based only on calculation shown in this paper. Hardy Cross method in
47
48 original form has historical value and should be replaced with Modified Hardy Cross method, or
49
50 even better with Node-Loop method.
51
52
53
54
55 Node-Loop method, presented among others in the text, is powerful numerical procedure for
56
57
58 calculation of flows in looped fluid distribution networks. Main advantage is that flow in each
59
60 pipe can be calculated directly, which is not possible after other available methods. In other
61
62 24
63
64
65
 methods, results of calculation are flow corrections which have to be added to flows calculated in
1
2 previous iteration using complex algebraic rules. Node-Loop method is recommended to be
3
4
5 used.
6
7
8
9
10 In the real network, consumers are not concentrated in the nodes. This can cause two-way flow
11
12 in some pipes (Brkić 2009). This can cause disturbance in convergence properties of certain
13
14
15
method. In such case, method should be changed. Some details on convergence properties can be
16
17 found in the paper of Mah (1974), Mah and Lin (1980) and Altman and Boulos (1995).
18
19 Simulation problem today can be solved using different software (Huddleston et al 2004,
20
21
22 Lopes 2004) in which can be implemented shown methods. Method for solution of pipe
23
24 equations proposed by Hamam and Brameller (1971) for gas networks and Todini and Pilati
25
26
27 (1988) for water networks is available, but not shown here.
28
29
30
31
32
Acknowledgements Financial support due scholarship program for PhD research given by the
33
34 Ministry of Science and Technological Development of Republic of Serbia is gratefully
35
36 acknowledged.
37
38
39
40
41 Notations:
42
43
44 p-pressure (Pa)
45
46 Δp-pressure drop (Pa)
47
48
49 C-pressure function in gas pipelines (Pa2)
50
51 F-pressure function C or Δp (Pa2 or Pa)
52
53
54
f- flow function (m3/s)
55
56 λ-Darcy friction factor (-)
57
58 L-pipe length (m)
59
60
61
62 25
63
64
65
 ρ-water density (kg/m3)
1
2 ρr-relative gas density (-)
3
4
5 Q-flow (m3/s)
6
7 D-pipe diameter (m)
8
9
10 Re-Reynolds number (-)
11
12 ε-pipe roughness (m)
13
14
15
ΔQ-correction of flow (m3/s)
16
17 Δx-correction of pressure function (Pa2 or Pa)
18
19 π~3.1415
20
21
22 subscripts
23
24 n-normal (~101325 Pa, ~274 K)
25
26
27 in-inner
28
29 r-relative
30
31
32
w-water
33
34 g-gas
35
36
37
38
39 References:
40
41 Altman T, Boulos PF (1995) Convergence of Newton method nonlinear network analysis.
42
43
44 Math Comput Model 21(4):35–41 doi:10.1016/0895-7177(95)00004-L
45
46 Andrijašev MM. (1964) Hydraulics calculation of water distribution networks. Stroizdat,
47
48
49 Moscow (in Russian)
50
51 Aynsley RM (1997) A resistance approach to analysis of natural ventilation airflow networks.
52
53
54
J Wind Eng Ind Aerod 67-68:711-719 doi:10.1016/S0167-6105(97)00112-8
55
56 Boulos PF, Lansey KE, Karney BW (2006) Comprehensive water distribution systems
57
58 analysis handbook for engineers and planners. MWH Soft, Hardback
59
60
61
62 26
63
64
65
 Brkić D (2009) An Improvement of Hardy Cross method applied on looped spatial natural gas
1
2 distribution networks. Appl Energ 86(7-8):1290-1300 doi:10.1016/j.apenergy.2008.10.005
3
4
5 Coelho PM, Pinho C (2007) Considerations about equations for steady state flow in natural
6
7 gas pipelines. J Braz Soc Mech Sci Eng 29(3):262–273 doi:10.1590/S1678-
8
9
10 58782007000300005
11
12 Colebrook CF (1939) Turbulent flow in pipes with particular reference to the transition region
13
14
15
between the smooth and rough pipe laws. J Inst Civil Eng (London) 11(4):133-156
16
17 doi:10.1680/ijoti.1939.13150
18
19 Corfield G, Hunt BE, Ott RJ, Binder GP, Vandaveer FE (1974) Distribution design for
20
21
22 increased demand. In: Segeler CG (ed) Gas Engineers Handbook, Industrial Press, New York,
23
24 pp 63–83
25
26
27 Cross H (1936) Analysis of flow in networks of conduits or conductors. Engineering
28
29 Experimental Station 286(34):3–29
30
31
32
Ekinci Ö, Konak H (2009) An optimization strategy for water distribution networks. Water
33
34 Resour Manag 23(1):169-185 doi:10.1007/s11269-008-9270-8
35
36 Epp R, Fowler AG (1970) Efficient code for steady flows in networks. J Hydraul Div ASCE
37
38
39 96(1):43–56
40
41 Farshad F, Rieke H, Garber J (2001) New developments in surface roughness measurements,
42
43
44 characterization, and modeling fluid flow in pipe. J Petrol Sci Eng 29(2):139–150
45
46 doi:10.1016/S0920-4105(01)00096-1
47
48
49 Gay B, Middleton P (1971) The solution of pipe network problems. Chem Eng Sci 26(1):109-
50
51 123 doi:10.1016/0009-2509(71)86084-0
52
53
54
Haaland SE (1983) Simple and explicit formulas for friction factor in turbulent pipe flow, J.
55
56 Fluid. Eng. T. ASME 105(1): 89-90 doi:10.1115/1.3240948
57
58
59
60
61
62 27
63
64
65
 Hamam YM, Brameller A (1971) Hybrid method for the solution of piping networks. Proc.
1
2 IEE 118(11):1607-1612 doi:10.1049/piee.1971.0292
3
4
5 Huddleston DH, Alarcon VJ, Chen W (2004) Water distribution network analysis using
6
7 Excel. J Hydraul Eng ASCE 130(10):1033-1035 doi:10.1061/(ASCE)0733-9429(2004)
8
9
10 Kumar SM, Narasimhan S, Bhallamudi SM (2010) Parameter estimation in water distribution
11
12 networks. Water Resour Manag 24(6):1251–1272 doi:10.1007/s11269-009-9495-1
13
14
15
Latišenkov AM, Lobačev VG (1956) Hydraulics. Gosstroizdat, Moscow (in Russian)
16
17 Lopes AMG (2004) Implementation of the Hardy-Cross method for the solution of piping
18
19 networks. Comput Appl Eng Educ 12(2):117-125 doi:10.1002/cae.20006
20
21
22 Mah RSH (1974) Pipeline network calculations using sparse computation techniques. Chem
23
24 Eng Sci 29(7):1629-1638 doi:10.1016/0009-2509(74)87014-4
25
26
27 Mah RSH, Lin TD (1980) Comparison of modified Newton‟s methods. Comput Chem Eng
28
29 4(2):75-78 doi:10.1016/0098-1354(80)80018-4
30
31
32
Mah RSH, Shacham M (1978) Pipeline network design and synthesis. Adv in Chem Eng
33
34 10:125-209 doi:10.1016/S0065-2377(08)60133-7
35
36 Mathews EH, Köhler PAJ (1995) A numerical optimization procedure for complex pipe and
37
38
39 duct network design. Int J Num Method Heat Fluid Flow 5(5):445-457
40
41 doi:10.1108/EUM0000000004072
42
43
44 Pretorius JJ, Malan AG, Visser JA (2008) A flow network formulation for compressible and
45
46 incompressible flow. Int J Num Method Heat Fluid Flow 18(2):185-201
47
48
49 doi:10.1108/09615530810846338
50
51 Shamir U, Howard CDD (1968) Water distribution systems analysis. J Hydraul Div ASCE
52
53
54
94:219–234
55
56
57
58
59
60
61
62 28
63
64
65
 Sukharev MG, Karasevich AM, Samoilov RV, Tverskoi IV (2005) Investigation of the
1
2 hydraulic resistance in polyethylene pipelines. J Eng Phys Thermophys 78(2):350-359
3
4
5 doi:10.1007/s10891-005-0068-8
6
7 Todini E, Pilati S (1988) A gradient method for the analysis of pipe networks. In “Computer
8
9
10 Applications for Water Supply and Distribution”, Vol. 1, John Wiley and Sons, pp. 1–20.
11
12 Wang Y-J, Hartman HL (1967) Computer solution of three-dimensional mine ventilation
13
14
15
networks with multiple fans and natural ventilation. Int J Rock Mech Min Sci 4(2):129-154
16
17 doi:10.1016/0148-9062(67)90039-3
18
19 Wood DJ, Charles COA (1972) Hydraulic network analysis using linear theory. J Hydraul
20
21
22 Div ASCE 98(7):1157–1170
23
24 Wood DJ, Rayes AG (1981) Reliability of algorithms for pipe network analysis. J Hydraul
25
26
27 Div ASCE 107(10):1145–1161
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62 29
63
64
65
 List of figures:
1
2 Figure 1. Example of pipeline network
3
4
5 Figure 2. Example network initial parameters prepared for loop oriented calculation
6
7 Figure 3. Contours for method of M.M. Andrijašev calculation
8
9
10 Figure 4. Example of pipeline network with three loops adjusted for node oriented methods
11
12 Figure 5. Comparisons of convergence for presented methods (gas network)
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62 30
63
64
65
 List of tables:
1
2 Table 1. Hardy Cross calculation for example gas network
3
4
5 Table 2. Hardy Cross calculation for example water network
6
7 Table 3. Calculation after the modified Hardy Cross method for example gas network
8
9
10 Table 4. Calculation after the modified Hardy Cross method for example water network
11
12 Table 5. Calculation after Modified M.M. Andrijašev method for example gas network
13
14
15
Table 6. Calculation after Modified M.M. Andrijašev method for example water network
16
17 Table 7. Calculation after Modified Node method for example gas network with three loops
18
19 Table 8. Calculation after Modified Node method for example water network with three loops
20
21
22 Table 9. Final flows for network presented in this paper
23
24 Table 10: Velocities for water and gas for calculated flows from example network
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62 31
63
64
65
line figure 1 DB
Click here to download high resolution image
line figure 2 DB
Click here to download high resolution image
colour figure 3 DB
Click here to download high resolution image
colour figure 4 DB
Click here to download high resolution image
colour figure 5 DB
Click here to download high resolution image
table 1 DB
Click here to download table: Table 1 DB.doc

Table 1. The Hardy Cross calculation for example gas network

Iteration 1 Iteration 2
F’= F’=

Loop Pipe a
Q b
C= p12  p22 C  Q  c
ΔQ1=
F d
ΔQ2 e
Q1=Q b
C= p12  p22 C  Q  c
ΔQ1=
F d
ΔQ2 Q2=Q
Q F' Q F'
3 -19444·10-5 -1264933339 11839776055 +13268·10-5 - -6177·10-5 -156904917 4623293467 -1128·10-5 - -7305·10-5
+9722·10 
-5
I 4 +2778·10 -5
+20357137 1333799622 +13268·10 -5 -5
+25767·10 -5 +1173109335 8285863414 -1128·10 -5572·10 = +19067·10-5
-5

7 -30556·10-5 -2399620963 14293015047 +13268·10-5 +881·10-5‡ -16407·10-5 -773797561 8583648143 -1128·10-5 -4154·10-5  -21689·10-5
Σ FI=-3644197165 27466590725 FI=+242406855 21492805024
1 +27778·10-5 +1344982709 8812326713 -9722·10-5 - +18056·10-5 +614087396 6189913452 +5572·10-5 - +23628·10-5
-5 -5 -5 -5
2 -27778·10 -200615476 1314432601 -9722·10 - -37500·10 -346390930 1681162091 +5572·10 - -31927·10-5
-13268·10 
II 4 -2778·10 -5
-20357137 1333799622 -9722·10 -5 -5
-25767·10 -5
-1173109335 8285863414 +5572·10 -5
+1128·10 ‡-5
-19067·10-5
+881·10  -4154·10 
-5
5 +2778·10 -5
+1828425 119798452 -9722·10 -5 -5
-6063·10 -5f -7569671 227216622.4 +5572·10 -5 -4645·10-5
Σ FII=+1125838521 11580357388 FII=-912982541 16384155579
5 -2778·10-5 -1828425 119798452 -881·10-5 +9722·10-5‡ +6063·10-5f +7569671 227216622 +4154·10-5 -5572·10-5= +4645·10-5
-5 -5 -5 -5
6 +2778·10 +65604940 4298435730 -881·10 - +1897·10 +32767180 3143915857 +4154·10 - +6051·10-5
III
+1128·10 
-5 -5 -5 -5 -5 -5 -5
7 +30556·10 +2399620963 14293015047 -881·10 -13268·10 = +16407·10 +773797561 8583648143 +4154·10 +21689·10
8 -16667·10-5 -2096864105 22897756035 -881·10-5 - -17548·10-5 -2302927589 23885524961 +4154·10-5 - -13394·10-5
Σ FIII=+366533372 41609005264 FIII=-1488793175 35840305583
Pipe lengths and diameters are shown in figure 1 and initial flow patterna in figure 2; 6 iterations are enough to achieve results shown in table 9 (this is for ΣF→0)
b
this is F calculated using (1), calso using (9), dΔQ2 is ΔQ1 from adjacent loop, efinal calculated flow in the first iteration is used for the calculation in the second iteration, fopposite flow direction than in
previous iteration
table 2 DB
Click here to download table: Table 2 DB.doc

Table 2. The Hardy Cross calculation for example water network

Iteration 1 Iteration 2
F’= F’=

Loop Pipe a
Q b
Δp=p1-p2 p  Q  c
ΔQ1=
F d
ΔQ2 e
Q1=Q b
Δp=p1-p2 p  Q  c
ΔQ1=
F d
ΔQ2 Q2=Q
Q F' Q F'
3 -19444·10-5 -5279095 54299272 +12306·10-5 - -7138·10-5 -750378 21023936 -400·10-5 - -7538·10-5
I 4 +2778·10-5 +69146 4978573 +12306·10-5 +9349·10-5  +24433·10-5 +4596449 37625131 -400·10-5 -4573·10-5= +19459·10-5
7 -30556·10-5 -10718549 70157780 +12306·10-5 +1068·10-5‡ -17181·10-5 -3450640 40167998 -400·10-5 -3920·10-5  -21501·10-5
Σ FI=-15928498 129435625 FI=+395430 98817066
1 +27778·10-5 +5919850 42622924 -9349·10-5 - +18429·10-5 +2639603 28646495 +4573·10-5 - +23002·10-5
2 -27778·10-5 -816585 5879413 -9349·10-5 - -37127·10-5 -1443784 7777597 +4573·10-5 - -32553·10-5
II -69146 4978573 -12306·10- -24433·10-5 -4596449 37625131 +4573·10-5 +400·10-5‡ -19459·10-5
4 -2778·10-5 -9349·10-5 5

5 +2778·10 -5 +5967 429677 -9349·10-5 +1068·10-5  -5503·10-5f -21346 775842 +4573·10-5 -0.03920  -4850·10-5
Σ FII=+5040086 53910587 FII=-3421977 74825067
5 -2778·10-5 -5967 429677 -1068·10-5 +9349·10-5‡ +5503·10-5f +21346 775842 +3920·10-5 -4573·10-5= +4850·10-5
6 +2778·10-5 +232186 16717415 -1068·10-5 - +1709·10-5 +92718 10848205 +3920·10-5 - +5629·10-5
III
+400·10  +21501·10
-5 -5 -5 -5 -5 -5 -5
7 +30556·10 +10718549 70157780 -1068·10 -12306·10 = +17181·10 +3450640 40167998 +3920·10
-5 -5 -5 -5
8 -16667·10 -8874257 106491089 -1068·10 - -17735·10 -10028129 113088168 +3920·10 - -13815·10-5
Σ FIII=+2070510 193795963 FIII=-6463424 164880214
Pipe lengths and diameters are shown in figure 1 and initial flow patterna in figure 2; 6 iterations are enough to achieve results shown in table 9 (this is for ΣF→0)
b
this is F calculated using (3), calso using (9) or (10), dΔQ2 is ΔQ1 from adjacent loop, efinal calculated flow in the first iteration is used for the calculation in the second iteration, fopposite flow direction
than in previous iteration
table 3 DB
Click here to download table: Table 3 DB.doc

Table 3. Calculation after the modified Hardy Cross method for example gas network
Iteration 1 Iteration 2
F’= F’=

Loop Pipe a
Q b
C= p p
2 2 C  Q  c
ΔQ1 d
ΔQ2 e
Q1=Q b
C= p p
2 2 C  Q  c
ΔQ1 d
ΔQ2 Q2=Q
1 2 1 2
Q Q
3 -19444·10-5 -1264933339 11839776055 +15114·10-5 - -4331·10-5 -82226369 3455539181 +2078·10-5 - -2253·10-5
+7937·10 
-5 -5 -5 -5 -5 -5 -5
I 4 +2778·10 +20357137 1333799622 +15114·10 +25828·10 +1178151884 8301891430 +2078·10 -5572·10 = +22263·10
7 -30556·10-5 -2399620963 14293015047 +15114·10-5 -4288·10-5  -19730·10-5 -1082435403 9985057887 +2078·10-5 -4154·10-5‡ -17620·10-5
Σ FI=-3644197166 27466590725 FI=+13490112 21742488499
1 +27778·10-5 +1344982710 8812326713 -7937·10-5 - +19841·10-5 +729037716 6687432899 +5643·10-5 - +25484·10-5
2 -27778·10-5 -200615476 1314432601 -7937·10-5 - -35715·10-5 -316967670 1615251716 +5643·10-5 - -30072·10-5
II
-15114·10-5  +1128·10-5 
-5
4 -2778·10 -5
-20357137 1333799622 -7937·10 -5
-25828·10 -5
-1178151884 8301891430 +5643·10 -22263·10-5
-5 -5 -5 -5f -5 -5
5 +2778·10 +1828426 119798452 -7937·10 -4288·10 = -9447·10 -16966070 326858346.3 +5643·10 -4154·10 ‡ -3772·10-5
Σ FII=+1125838522 11580357388 FII=-783047908 16931434392
5 -2778·10-5 -1828425 119798452 +4288·10-5 +7937·10-5‡ +9447·10-5f +16966070 326858346.3 -32·10-5 -5572·10-5= +3772·10-5
6 +2778·10-5 +65604940 4298435730 +4288·10-5 - +7066·10-5 +358812846 9242405949 -32·10-5 - +7034·10-5
III -5 -5 -5 -5 -5
7 +30556·10 +2399620963 14293015047 +4288·10 -15114·10 = +19730·10 +1082435403 9985057887 -32·10 +1128·10 = +17620·10-5
-5
-5 -5 -5 -5
8 -16667·10 -2096864105 22897756035 +4288·10 - -12379·10 -1220331619 17942054130 -32·10 - -12411·10-5
Σ FIII=+366533372 41609005264 FIII=+237882701 37496376312
Pipe lengths and diameters are shown in figure 1 and initial flow patterna in figure 2; 3 iterations are enough to achieve results shown in table 9 (this is for ΣF→0)
b
this is F calculated using (1), cusing (12) where (Δp) is replaced with (C), dΔQ2 is ΔQ1 from adjacent loop, efinal calculated flow in the first iteration is used for the calculation in the second iteration, fopposite
flow direction than in previous iteration
table 4 DB
Click here to download table: Table 4 DB.doc

Table 4. Calculation after the modified Hardy Cross method for example water network
Iteration 1 Iteration 2
F’= F’=

Loop Pipe a
Q b
Δp=p1-p2 p  Q  c
ΔQ1 d
ΔQ2 e
Q1=Q b
Δp=p1-p2 p  Q  c
ΔQ1 d
ΔQ2 Q2=Q
Q Q
3 -19444·10-5 -5279095 54299272 +14195·10-5 - -5250·10-5 -416024 15849700 +2505·10-5 - -2744·10-5
+8006·10 
-5 -5 -5
I 4 +2778·10 -5 +69146 4978573 +14195·10 -5 -5
+24978·10 +4800902 38440414 +2505·10 -5312·10 = +22172·10-5
-4053·10  -107·10 
-5 -5
7 -30556·10 -5
-10718549 70157780 +14195·10 -5 -5 -20413·10 -4840373 47423556 +2505·10 -5 -18015·10-5
Σ FI=-15928498 129435625 FI=-455494 101713670
1 +27778·10-5 +5919850 42622924 -8006·10-5 - +19772·10-5 +3030676 30656275 +5312·10-5 - +25084·10-5
-5 -5 -5 -5
2 -27778·10 -816585 5879413 -8006·10 - -35784·10 -1342790 7505063 +5312·10 - -30472·10-5
-14195·10  -2505·10 
II -5
-69146 4978573 -5 -5 -24978·10 -5
-4800902 38440414 +5312·10 -5 -5 -22172·10-5
4 -2778·10 -8006·10
-107·10-5 
-5 -5
5 +2778·10-5 +5967 429677 -8006·10-5 -4053·10-5= -9281·10 -57430 1237632 +5312·10 -4076·10-5
Σ FII=+5040086 53910587 FII=-3170446 77839384
5 -2778·10-5 -5967 429677 +4053·10-5 +8006·10-5‡ +9281·10-5 +57430 1237632 +107·10-5 -5312·10-5= +4076·10-5
-5 -5 -5 -5
6 +2778·10 +232186 16717415 +4053·10 - +6830·10 +1312311 38425447 +107·10 - +6937·10-5
III -5 -5 -5 -5 -5 -5
7 +30556·10 +10718549 70157780 +4053·10 -14195·10 = +20413·10 +4840373 47423556 +107·10 -2505·10 = +18015·10-5
8 -16667·10-5 -8874257 106491089 +4053·10-5 - -12614·10-5 -5136681 81444010 +107·10-5 - -12507·10-5
Σ FIII=+2070510 193795963 FIII=+1073434 168530645
Pipe lengths and diameters are shown in figure 1 and initial flow patterna in figure 2; 3 iterations are enough to achieve results shown in table 9 (this is for ΣF→0)
b
this is F calculated using (3), cusing (12), dΔQ2 is ΔQ1 from adjacent loop, efinal calculated flow in the first iteration is used for the calculation in the second iteration, fopposite flow direction than in
previous iteration
table 5 DB
Click here to download table: Table 5 DB.doc

Table 5. Calculation after Modified M.M. Andrijašev method for example gas network
Iteration 1 Iteration 2
F’= F’=

ContourO Pipe a
Q b
C= p12  p22 C  Q  c
Q1o d
Q2o e
Q1=Q b
C= p12  p22 C  Q  c
Q1o d
Q2o Q2=Q
Q Q
6 +2778·10-5 65604941 4298435730 +13669·10-5 -9381·10-5‡ +7066·10-5 358812846 9242405949 -1798·10-5 +1767·10-5  +7034·10-5
8 -16667·10-5 -2096864106 22897756035 +13669·10-5 -9381·10-5  -12379·10-5 -1220331619 17942054130 -1798·10-5 +1767·10-5= -12411·10-5
IO 3 -19444·10-5 -1264933339 11839776055 +13669·10-5 +1444·10-5= -4331·10-5 -82226369 3455539181 -1798·10-5 +3876·10-5= -2253·10-5
4 +2778·10-5 20357137 1333799622 +13669·10-5 +9381·10-5  +25828·10-5 1178151884 8301891430 -1798·10-5 -1767·10-5= +22263·10-5
5 -2778·10-5 -1828426 119798452 +13669·10-5 -1444·10-5  +9447·10-5f 16966070 326858346 -1798·10-5 -3876·10-5= +3772·10-5
Σ F°=-3277663793 40489565894 F°=251372813 39268749037
6 +2778·10-5 65604941 4298435730 -9381·10-5 +13669·10-5  +7066·10-5 358812846 9242405949 +1767·10-5 -1798·10-5‡ +7034·10-5
8 -16667·10-5 -2096864106 22897756035 -9381·10-5 +13669·10-5= -12379·10-5 -1220331619 17942054130 +1767·10-5 -1798·10-5  -12411·10-5
7 +30556·10-5 2399620963 14293015047 -9381·10-5 -1444·10-5= +19730·10-5 1082435403 9985057887 +1767·10-5 -3876·10-5= +17620·10-5
IIO
4 -2778·10-5 -20357137 1333799622 -9381·10-5 -13669·10-5  -25828·10-5 -1178151884 8301891430 +1767·10-5 +1798·10-5‡ -22263·10-5
2 -27778·10-5 -200615476 1314432601 -9381·10-5 +1444·10-5= -35715·10-5 -316967670 1615251716 +1767·10-5 +3876·10-5= -30072·10-5
1 +27778·10-5 1344982710 8812326713 -9381·10-5 +1444·10-5  +19841·10-5 729037716 6687432899 +1767·10-5 +3876·10-5  +25484·10-5
Σ F°=1492371895 52949765748 F°=-545165207 53774094012
1 +27778·10-5 1344982710 8812326713 +1444·10-5 -9381·10-5‡ +19841·10-5 729037716 6687432899 +3876·10-5 +1767·10-5  +25484·10-5
5 +2778·10-5 1828426 119798452 +1444·10-5 -13669·10-5= -9447·10-5f -16966070 326858346 +3876·10-5 +1798·10-5‡ -3772·10-5
IIIO 7 -30556·10-5 -2399620963 14293015047 +1444·10-5 +9381·10-5‡ -19730·10-5 -1082435403 9985057887 +3876·10-5 -1767·10-5  -17620·10-5
3 -19444·10-5 -1264933339 11839776055 +1444·10-5 +13669·10-5= -4331·10-5 -82226369 3455539181 +3876·10-5 -1798·10-5  -2253·10-5
2 -27778·10-5 -200615476 1314432601 +1444·10-5 -9381·10-5  -35715·10-5 -316967670 1615251716 +3876·10-5 +1767·10-5= -30072·10-5
Σ F°=-2518358644 36379348868 F°=-769557797 22070140031
Pipe lengths and diameters are shown in figure 1, initial flow pattern in figure 3 and numerical values for initial flowsa in figure 2; 3 iterations are enough to achieve results shown in table 9 (this is for ΣF→0)
b
this is F calculated using (1), cusing (15),
d
Q2o is Q1o from adjacent contour, efinal calculated flow in the first iteration is used for the calculation in the second iteration, fopposite flow direction than in
previous iteration
table 6 DB
Click here to download table: Table 6 DB.doc

Table 6. Calculation after Modified M.M. Andrijašev method for example water network
Iteration 1 Iteration 2
F’= F’=

ContourO Pipe a
Q b
Δp=p1-p2 p  Q  c
Q1o d
Q2o e
Q1=Q b
Δp=p1-p2 p  Q  c
Q1o d
Q2o Q2=Q
Q Q
6 +2778·10-5 +232186 16717415 +13127·10-5 -9074·10-5‡ +6830·10-5 +1312311 19212723 -1350·10-5 +1457·10-5  +6937·10-5
8 -16667·10-5 -8874257 106491089 +13127·10-5 -9074·10-5  -12614·10-5 -5136681 40722005 -1350·10-5 +1457·10-5= -12507·10-5
IO 3 -19444·10-5 -5279095 54299272 +13127·10-5 +1068·10-5= -5250·10-5 -416024 7924850 -1350·10-5 +3855·10-5= -2744·10-5
4 +2778·10-5 +69146 4978573 +13127·10-5 +9074·10-5  +24978·10-5 4800902 19220207 -1350·10-5 -1457·10-5= +22172·10-5
5 -2778·10-5 -5967 429677 +13127·10-5 -1068·10-5  +9281·10-5f +57430 618816 -1350·10-5 -3855·10-5= +4076·10-5
Σ F°=-13857988 182916028 F°=617939 87698601
6 +2778·10-5 +232186 16717415 -9074·10-5 +13127·10-5  +6830·10-5 1312311 19212723 +1457·10-5 -1350·10-5‡ +6937·10-5
8 -16667·10-5 -8874257 106491089 -9074·10-5 +13127·10-5= -12614·10-5 -5136681 40722005 +1457·10-5 -1350·10-5  -12507·10-5
7 +30556·10-5 +10718549 70157780 -9074·10-5 -1068·10-5= +20413·10-5 4840373 23711778 +1457·10-5 -3855·10-5= +18015·10-5
IIO
4 -2778·10-5 -69146 4978573 -9074·10-5 -13127·10-5  -24978·10-5 -4800902 19220207 +1457·10-5 +1350·10-5‡ -22172·10-5
2 -27778·10-5 -816585 5879413 -9074·10-5 +1068·10-5= -35784·10-5 -1342790 3752532 +1457·10-5 +3855·10-5= -30472·10-5
1 +27778·10-5 +5919850 42622924 -9074·10-5 +1068·10-5  +19772·10-5 3030676 15328137 +1457·10-5 +3855·10-5  +25084·10-5
Σ F°=+7110597 246847196 F°=-2097013 121947382
1 +27778·10-5 +5919850 42622924 +1068·10-5 -9074·10-5‡ +19772·10-5 3030676 15328137 +3855·10-5 +1457·10-5  +25084·10-5
5 +2778·10-5 +5967 429677 +1068·10-5 -13127·10-5= -9281·10-5f -57430 618816 +3855·10-5 +1350·10-5‡ -4076·10-5
IIIO 7 -30556·10-5 -10718549 70157780 +1068·10-5 +9074·10-5‡ -20413·10-5 -4840373 23711778 +3855·10-5 -1457·10-5  -18015·10-5
3 -19444·10-5 -5279095 54299272 +1068·10-5 +13127·10-5= -5250·10-5 -416024 7924850 +3855·10-5 -1350·10-5  -2744·10-5
2 -27778·10-5 -816585 5879413 +1068·10-5 -9074·10-5 
-35784·10-5 -1342790 3752532 +3855·10-5 +1457·10-5=
-30472·10-5
Σ F°=-10888413 173389067 F°=-3625941 51336113
Pipe lengths and diameters are shown in figure 1, initial flow pattern in figure 3 and numerical values for initial flowsa in figure 2; 3 iterations are enough to achieve results shown in table 9 (this is for ΣF→0)
b
this is F calculated using (3), cusing (15),
d
Q2o is Q1o from adjacent contour, efinal calculated flow in the first iteration is used for the calculation in the second iteration, fopposite flow direction than in
previous iteration
table 7 DB
Click here to download table: Table 7 DB.doc

Table 7. Calculation after Modified Node method for example gas network with three loops
Iteration 1 Iteration 2
a b c d
Node Pipe C (Pa2) a
Q=f(C) b
f’ c
ΔC1 d
ΔC2 e
C1=C Q=f(C) f’ ΔC1 ΔC2 C2=C
2 +0.25·1010 +1.1108 24415·10-14 -36.1·108 - -11.1·108f -0.7109 35204·10 -14
+28.4·10 8
- +17.3·108f
I 3 +0.25·1010 +0.2827 6214·10-14 -36.1·108 -6.99·108= -18.1·108f -0.2367 7189·10-14 +28.4·108 +4.57·108‡ +14.9·108f
4 -0.50·1010 -0.5715 6280·10-14 -36.1·108 +97.0·108‡ +10.9·108f +0.2475 12474·10-14 +28.4·108 -66.9·108‡ -27.6·108f
Constant output flow -0.0555 -0.0555
Σf=+0.7665 36909·10-14 Σf=-0.7557 54868·10-14
3 -0.25·1010 -0.2827 6214·10-14 +6.99·108 +36.1·108‡ +18.1·108f +0.2367 7189·10-14 -4.57·108 -28.4·108= -14.9·108f
II 7 -0.75·1010 -0.5715 4187·10-14 +6.99·108 +97.0·108‡ +29.0·108f +0.3390 6425·10-14 -4.57·108 -66.9·108= -42.5·108f
8 -0.75·1010 -0.3357 2459·10-14 +6.99·108 +78.0·108‡ +9.99·108f +0.1109 6099·10-14 -4.57·108 -39.4·108= -34.0·108f
Constant input flow +0.2777 +0.2777
Σf=-0.9121 12860·10-14 Σf=+0.9644 19713·10-14
4 +0.50·10 10
+0.5715 6280·10-14 -97.0·10 8 +36.1·10 8
 -10.9·10 8f -0.2475 12474·10-14 +66.9·108 -28.4·108  +27.6·108f
III 5 +0.50·1010 +2.1483 23608·10-14 -97.0·108 +7.52·10 8
 -39.5·108f -1.8867 26260·10
-14
+66.9·108 5.29·108‡ +32.7·108f
7 +0.75·1010 +0.5715 4187·10-14 -97.0·108 -6.99·108= -29.0·108f -0.3390 6425·10 -14
+66.9·10 8
4.57·108‡ +42.5·108f
Constant output flow -0.3611 -0.3611
Σf=+2.9302 34075·10-14 Σf=-2.8343 45160·10-14
10 -14 8 8
1 +0.25·10 +0.3905 8582·10 -7.52·10 - +17.5·10 +0.3208 10084·10-14 -5.29·108 - +12.2·108
IV 5 -0.50·1010 -2.1483 23608·10-14 -7.52·108 +97.0·108‡ +39.5·108f +1.8867 26260·10-14 -5.29·108 -66.9·108= -32.7·108f
10 -14 8 8 8f -14 8 8
6 -0.50·10 -0.3004 3302·10 -7.52·10 +78.0·10 ‡ +20.5·10 +0.1840 4937·10 -5.29·10 -39.4·10 = -24.3·108f
Constant output flow -0.2222 -0.2222
Σf=-2.2805 35492·10-14 Σf=2.1693 41281·10-14
+7.52·108 
-14
6 +0.50·1010 +0.3004 3302·10-14 -78.0·108 -20.5·108f -0.1840 4937·10 39.4·108 5.29·108‡ +24.3·108f
V 10 -14 8 8 8f -14 8 8
8 +0.75·10 +0.3357 2459·10 -78.0·10 -6.99·10 = -9.99·10 -0.1109 6099·10 39.4·10 4.57·10 ‡ +34.0·108f
Constant output flow -0.1944 -0.1944
Σf=+0.4417 5761·10-14 Σf=-0.4893 11036·10-14
Pipe lengths and diameters are shown in figure 1; see figure 4 for initial pattern (red letters); 26 iterations are enough to achieve results shown in table 9 (this is for Σf→0)
a
using (17), bF’=|∂Q(C)/∂C|,cΔC1 after eq. (20), dΔC2 is ΔC1 from adjacent node, efinal calculated pressure function in the first iteration is used for the calculation in the second iteration, fopposite
flow direction than in previous iteration
table 8 DB
Click here to download table: Table 8 DB.doc

Table 8. Calculation after Modified Node method for example water network with three loops
Iteration 1 Iteration 2
a b c d
Node Pipe Δp (Pa) a
Q=f(Δp) b
f’ c
ΔΔp1 d
ΔΔp2 e
Δp1=Δp Q=f(Δp) f’ ΔΔp1 ΔΔp2 Δp2=Δp
2 +25·103 +4953·10 -5
9907·10 -10
+177553 - +202553 +12974·10 -5
3203·10 -10
+552198 - +754750
I 3 +25·103 +1360·10-5 2720·10-10 +177553 -167854= +34699 +1396·10-5 2010·10-10 +552198 -730668= -143771f
-798781  -2779111 
-5 -10
4 -50·103 -2580·10-5 2500·10-10 +177553 -671229 -8616·10 642·10 +552198 -2898142
Constant output flow -5556·10-5 -5556·10-5
Σf -1823·10-5 15207·10-10 +198·10-5 5856·10-10
-177553  -552198 
-5
3 -25·10 3
-1360·10 -5
2720·10 -10
+167854 -34699 -1396·10 2010·10-10 +730668 +143771f
-798781  -2779111 
3 -5 -10 -5 -10
II 7 -75·10 -2580·10 1720·10 +167854 -705927 -7214·10 511·10 +730668 -2754371
8 -75·103 -1555·10-5 1037·10-10 +167854 -1244652  -1151798 -5453·10-5 237·10-10 +730668 -3619538  -4040669
Constant input flow 27778·10-5 +27778·10-5
Σf 22282·10-5 5477·10-10 +13715·10-5 2759·10-10
4 +50·103 +2580·10-5 2580·10-10 +798781 -177553= +671229 +8616·10-5 642·10-10 +2779111 -552198= +2898142
-2785111 
3 -5 -10 f -5 -10
III 5 +50·10 +9206·10 9206·10 +798781 -854621= -5840 -2958·10 25329·10 +2779111 -11840
7 +75·103 +2580·10-5 1720·10-10 +798781 -167854= +705927 +7214·10-5 511·10-10 +2779111 -730668= +2754371
Constant output flow -36111·10-5 -36111·10-5
Σf -21744·10-5 13507·10-10 -23239·10-5 26482·10-10
1 +25·103 +1825·10-5 3649·10-10 +854621 - +879621 +9639·10-5 548·10-10 +2785111 - +3664732
-798781 
3 -5 -10 f -5 -10
IV 5 -50·10 -9206·10 9206·10 +854621 +5840 +2958·10 25329·10 +2785111 -2779111= +11840
-1244652  -3619538 
3 -5 -10 -5 -10
6 -50·10 -1372·10 1372·10 +854621 -440031 -3672·10 417·10 +2785111 -1274458
Constant output flow -22222·10-5 -22222·10-5
Σf -30976·10-5 14228·10-10 -13297·10-5 26294·10-10
3 -5 -10 -5
6 +50·10 +1372·10 1372·10 +1244652 -854621= +440031 +3672·10 417·10-10 +3619538 -2785111= +1274458
V 3 -5 -10 -5
8 +75·10 +1555·10 1037·10 +1244652 -167854= +1151798 +5453·10 237·10-10 +3619538 -730668= +4040669
Constant output flow -19444·10-5 -19444·10-5
Σf -16517·10-5 2409·10-10 -10320·10-5 654·10-10
Pipe lengths and diameters are shown in figure 1; see figure 4 for initial pattern (blue letters); 9 iterations are enough to achieve results shown in table 9 (this is for Σf→0)
a
using (18); to calculate friction factor λ, Reynolds number has to be calculated and for this velocity have to be chosen (this velocity does not have effect on final results, here is chosen extremely
large velocity 100 m/s), bF’=|∂Q(Δp)/∂Δp|,cΔΔp1 after eq. (20), dΔΔp2 is ΔΔp1 from adjacent node, efinal calculated pressure function in the first iteration is used for the calculation in the second
iteration, fopposite flow direction than in previous iteration
table 9 DB
Click here to download table: Table 9 DB.doc

Table 9. Final flows for network presented in this paper

Flows (m3/h)
Pipe 1 2 3 4 5 6 7 8
Water 902.27 1097.73 94.86 802.87 -146.23a 248.50 643.36 451.50
Gas 913.72 1086.28 82.01 804.27 -137.86a 251.58 633.60 448.42
a
sing minus means flow direction opposite than first assumed in figure 2
table 10 DB
Click here to download table: Table 10 DB.doc

Table 10: Velocities for water and gas for calculated flows from example network
Velocity (m/s)
Pipe 1 2 3 4 5 6 7 8
Water 21.0 15.4 2.2 18.6 -1.7a 9.4 14.9 13.2
Gas 5.3 3.8 0.5 4.7 -0.4a 2.4 3.7 3.3
a
sing minus means flow direction opposite than first assumed in figure 2