"""
Demostration of the Hardy-Cross relaxation solution of the pipeline network
(1) (2) (3) (7)
----o-----o-----o------
\(4) |(5) /(6)
\ | /
\--o--/
"""
import numpy as n
import pylab as p
L = n.array([2,2,3,3,2,3,2],dtype='f') # km
D = n.array([25,25,20,20,15,20,25],dtype='f') # cm
C = n.array([120,100,110,120,130,120,130],dtype='f') # Hazen Williams
# resistance per unit length as a function
r = lambda l,d,c: 1.526e7/(c**1.852*d**4.87)*l
# head loss as a function:
hf = lambda R,Q: R*n.sign(Q)*n.abs(Q)**1.852
R = r(L,D,C)# resistances
# define branches - each row contains numbers of pipes:
branch = n.array([[2,5,4],[3,6,5]])-1 # Python counts from zero, not 1
rows,cols = branch.shape # rows = num of branches, cols = pipes in each
# initial guess that
Q = n.array([500,250,250,-250,0,-250,500],dtype='f') # m**3/hr
dQ = 500.0
# main loop
while abs(dQ) > 0.5:
for i in n.arange(rows):
y = hf(R,Q)
yq = n.abs(1.852*y/n.abs(Q))
yq[n.isnan(yq)] = 0.0
# for the first branch
sumyq = n.sum(yq[branch[i,:]])
sumy = n.sum(y[branch[i,:]])
dQ = -1*sumy/sumyq
print("dQ = %f" % dQ)
Q[branch[i,:]] += dQ
print("Discharges [m^3/hr]:\n")
print Q
# we're looking for equivalent pipe solution with:
Deq = 25 # cm
Ceq = 120 #
y = hf(R[1],Q[1])+hf(R[2],Q[2])
Leq = y/(r(1,25,120)*Q[0]**1.852)
Leq = Leq + L[0] + L[6]*(120/C[6])**1.852
print("Equivalent length = %3.2f km" % Leq)
|
|
Warning: invalid value encountered in divide
dQ = 70.284711
dQ = -18.083459
dQ = -6.477074
dQ = 0.765176
dQ = -0.338565
dQ = 0.067968
Discharges [m^3/hr]:
[ 500. 313.46908569 232.74967957 -186.53092957 46.21875763
-267.25033569 500. ]
Equivalent length = 7.44 km
Warning: invalid value encountered in divide
dQ = 70.284711
dQ = -18.083459
dQ = -6.477074
dQ = 0.765176
dQ = -0.338565
dQ = 0.067968
Discharges [m^3/hr]:
[ 500. 313.46908569 232.74967957 -186.53092957 46.21875763
-267.25033569 500. ]
Equivalent length = 7.44 km
|