流体网络拓扑（3）——网络分流_流体网络理论

5网络分流

5.1网络分流预处理

流体网络的源点与汇点之间保持流量平衡，从源点流入网络的总流量等于从汇点流出网络的总流量，只是源点与汇点的压能不同。在多源多汇的情况下，源点之间以及汇点之间的压能有时也是不一样的。为了使求解模型具有通用性，本书通过增加虚拟节点和虚拟分支的办法将有源汇的流体网络变成无源汇的流体网络。
视具体情况，对单一源汇网络可以用一虚拟分支直接将源汇节点连接起来，也可以把源点和汇点直接连接起来。对多源多汇网络也可以直接将源汇节点合并。

5.2流体网络的基本定律

5.2.1质量守恒定律

5.2.2能量守恒定律

5.2.3阻力定律

例1

5.3网络分流算法

5.3.1网络分流算法综述

当流体网络中所有的流阻为已知，并已知网络的总流量、或已知回路的附加阻力，求所有分支流量的过程叫做网络分流，也称网络解算。
如果将分支流量视为未知数，那么通过基本关联矩阵可以构造m-1个线性独立的关于分支流量的方程组；如果把分支的流阻作为已知数，那么利用基本回路矩阵也可以构造n-m+1个关于分支流量的二次非线性方程组。将流量平衡方程和阻力平衡方程联立，n个分支，（m-1）+（n-m+1）= n 个方程，一般来说理论上是有唯一解的。
流体网络分流算法研究最早可以追溯到1854年，J.J.Atkinson在北英格兰采矿工程师学会上发表的一篇论文奠定了一般流体网络分流理论的基础。1928年波兰学者H.Czecaott提出了通风网络解算问题，1925年发表了关于通风网络中的角联分支理论，1931年提出了用几何学的方法解含有角联分支的θ型通风网络。1935年波兰学者S.Barczyk根据非线性方程组的Newton解法，提出了一种通风网络的近似解算方法。1936年美国学者H.Cross提出一种解算流体管道网络的逐次计算法。1951年英国学者D.Scott、F.Hinsley在英国的《煤炭工程》上发表了“通风网络理论”一文，文中系统提出了通风网络的解算方法，该法在本质上与Cross法相同，也属于逐次迭代计算方法。此外，1938年英国学者S.Weeks尝试了用图解法来解算简单的通风网络，1950年美国学者M.Mellrog、荷兰学者W.Mass采用钨丝灯泡来模拟通风网络，日本学者在这段时间相继提出了京大一法和京大二法求解通风网络。
综上，网络解法可分为：解析法、图解法、物理相似模拟法、数值方法。数值法属于近似法，是目前研究分流的主要手段。从计算数学的角度看，数值方法可分为3类：斜量法、迭代法和直接代入法。本节主要介绍Barczyk法和Cross法，它们分别属于斜量法和迭代法。

5.3.2 Barczyk法

例2

无回路附加阻力的Barczyk法迭代过程程序

#Barczyk
from sympy import *

q1,r1 = 20,0.375
q2,r2 = 11.5,1.5
q3,r3 = 8.5,2
q4,r4 = 10.5,4.6875
q5,r5 = 1,12.5
q6,r6 = 9.5,1.7361
q7,r7 = 20,0.5
q8 = 20

C_T = Matrix([[0,-1,-1,0,0],[0,0,-1,-1,0],[1,1,1,1,1]])
i = 0
while True:
    i += 1
    print(f'第{
     i}次迭代：')
    f1 = r3 * q3 * abs(q3) - r2 * q2 * abs(q2) - r5 * q5 * abs(q5)
    f2 = r4 * q4 * abs(q4) - r5 * q5 * abs(q5) - r6 * q6 * abs(q6)

    df13 = 2*r3*abs(q3)+2*r2*abs(q2)+2*r5*abs(q5)
    df14 = 0
    df23 = 0
    df24 = 2*r4*abs(q4)+2*r5*abs(q5)+2*r6*abs(q6)
    print('f1=',f1,'f2=',f2)
    print('df13=',df13,'df24=',df24)
    J = Matrix([[df13,df14],[df23,df24]])
    print('J:',J)
    J_inv = J**(-1)
    print('J的逆矩阵:',J_inv)
    F = Matrix([[f1],[f2]])
    print('F：',F)
    var_Q = -1*J_inv*F
    print('Q的变量',var_Q)
    max_vq = max(abs(var_Q[0]),abs(var_Q[1]))
    max_f = max(abs(f1),abs(f2))
    if max_vq <= 10**-3 and max_f <= 10**-3:
        print(q1,q2,q3,q4,q5,q6)
        print(f'max△q = {
     max_vq} , max△f = {
     max_f}')
        break
    else:
        q3 = q3 + var_Q[0]
        q4 = q4 + var_Q[1]
        Q = Matrix([[q3,q4,q8]]) * C_T
        q1,q2,q5,q6,q7 = Q[0],Q[1],Q[2],Q[3],Q[4]
        print(q1, q2, q3, q4, q5, q6)
        print(f'max△q = {
     max_vq} , max△f = {
     max_f}')
1
2
3
4
5
6
7
8
9
10
11
12
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

结果：

第1次迭代：
f1= -66.375 f2= 347.61384999999996
df13= 93.5 df24= 156.42340000000002
J: Matrix([[93.5000000000000, 0], [0, 156.423400000000]])
J的逆矩阵: Matrix([[0.0106951871657754, 0], [0, 0.00639290540929298]])
F： Matrix([[-66.3750000000000], [347.613850000000]])
Q的变量 Matrix([[0.709893048128342], [-2.22226246201016]])
20 10.7901069518717 9.20989304812834 8.27773753798984 2.51236941388182 11.7222624620102
max△q = 2.22226246201016 , max△f = 347.61384999999996
第2次迭代：
f1= -83.8953530309290 f2= 3.73190331094872
df13= 132.019128395174 df24= 181.115064486292
J: Matrix([[132.019128395174, 0], [0, 181.115064486292]])
J的逆矩阵: Matrix([[0.00757465991599864, 0], [0, 0.00552135187007422]])
F： Matrix([[-83.8953530309290], [3.73190331094872]])
Q的变量 Matrix([[0.635478767741932], [-0.0206051513248429]])
20 10.1546281841297 9.84537181587027 8.25713238666500 1.89749579746473 11.7428676133350
max△q = 0.635478767741932 , max△f = 83.8953530309290
第3次迭代：
f1= -5.81814671887107 f2= 35.1893190116522
df13= 117.282766752488 df24= 165.621595988624
J: Matrix([[117.282766752488, 0], [0, 165.621595988624]])
J的逆矩阵: Matrix([[0.00852640185501748, 0], [0, 0.00603785994230296]])
F： Matrix([[-5.81814671887107], [35.1893190116522]])
Q的变量 Matrix([[0.0496078569765461], [-0.212468179657375]])
20 10.1050203271532 9.89497967284682 8.04466420700762 2.06035612014556 11.9553357929924
max△q = 0.212468179657375 , max△f = 35.1893190116522
第4次迭代：
f1= -10.4092500389306 f2= 2.15495819011008
df13= 121.403882676486 df24= 168.438946884764
J: Matrix([[121.403882676486, 0], [0, 168.438946884764]])
J的逆矩阵: Matrix([[0.00823696885102742, 0], [0, 0.00593686922469388]])
F： Matrix([[-10.4092500389306], [2.15495819011008]])
Q的变量 Matrix([[0.0857406683332272], [-0.0127937049593666]])
20 10.0192796588200 9.98072034118005 8.03187050204826 1.98740915677170 11.9681294979517
max△q = 0.0857406683332272 , max△f = 10.4092500389306
第5次迭代：
f1= -0.721829720024552 f2= 4.35037510719090
df13= 119.665949260472 df24= 166.539754118783
J: Matrix([[119.665949260472, 0], [0, 166.539754118783]])
J的逆矩阵: Matrix([[0.00835659605911233, 0], [0, 0.00600457233344274]])
F： Matrix([[-0.721829720024552], [4.35037510719090]])
Q的变量 Matrix([[0.00603203939370732], [-0.0261221420087365]])
20 10.0132476194262 9.98675238057375 8.00574836003952 2.00749925938673 11.9942516399605
max△q = 0.0261221420087365 , max△f = 4.35037510719090
第6次迭代：
f1= -1.30291156560526 f2= 0.296672041245159
df13= 120.174233865242 df24= 166.887812904309
J: Matrix([[120.174233865242, 0], [0, 166.887812904309]])
J的逆矩阵: Matrix([[0.00832125130185024, 0], [0, 0.00599204928506902]])
F： Matrix([[-1.30291156560526], [0.296672041245159]])
Q的变量 Matrix([[0.0108418545614885], [-0.00177767349264302]])
20 10.0024057648648 9.99759423513524 8.00397068654688 1.99843507831788 11.9960293134531
max△q = 0.0108418545614885 , max△f = 1.30291156560526
第7次迭代：
f1= -0.0901851748231124 f2= 0.543107709620699
df13= 119.958471193082 df24= 166.650715126496
J: Matrix([[119.958471193082, 0], [0, 166.650715126496]])
J的逆矩阵: Matrix([[0.00833621827666030, 0], [0, 0.00600057431041296]])
F： Matrix([[-0.0901851748231124], [0.543107709620699]])
Q的变量 Matrix([[0.000751803302644234], [-0.00325895817013719]])
20 10.0016539615621 9.99834603843789 8.00071172837674 2.00094223318537 11.9992882716233
max△q = 0.00325895817013719 , max△f = 0.543107709620699
第8次迭代：
f1= -0.162898698364152 f2= 0.0375135257347381
df13= 120.021901868072 df24= 166.694157019896
J: Matrix([[120.021901868072, 0], [0, 166.694157019896]])
J的逆矩阵: Matrix([[0.00833181264782155, 0], [0, 0.00599901051049222]])
F： Matrix([[-0.162898698364152], [0.0375135257347381]])
Q的变量 Matrix([[0.00135724143534411], [-0.000225044035168314]])
20 10.0002967201268 9.99970327987323 8.00048668434157 1.99981003578520 11.9995133156584
max△q = 0.00135724143534411 , max△f = 0.162898698364152
第9次迭代：
f1= -0.0112726051922394 f2= 0.0678781688012862
df13= 119.994954174503 df24= 166.664523694961
J: Matrix([[119.994954174503, 0], [0, 166.664523694961]])
J的逆矩阵: Matrix([[0.00833368375261635, 0], [0, 0.00600007714797335]])
F： Matrix([[-0.0112726051922394], [0.0678781688012862]])
Q的变量 Matrix([[9.39423267402240e-5], [-0.000407274249470875]])
20 10.0002027778000 9.99979722219997 8.00007941009210 2.00012336770793 11.9999205899079
max△q = 0.000407274249470875 , max△f = 0.0678781688012862
第10次迭代：
f1= -0.0203630010839220 f2= 0.00469593253927769
df13= 120.002881414898 df24= 166.669952934590
J: Matrix([[120.002881414898, 0], [0, 166.669952934590]])
J的逆矩阵: Matrix([[0.00833313323988112, 0], [0, 0.00599988169668742]])
F： Matrix([[-0.0203630010839220], [0.00469593253927769]])
Q的变量 Matrix([[0.000169687601196165], [-2.81750396912911e-5]])
20 10.0000330901988 9.99996690980117 8.00005123505241 1.99998185514642 11.9999487649476
max△q = 0.000169687601196165 , max△f = 0.0203630010839220
第11次迭代：
f1= -0.00140907480732722 f2= 0.00848465542944155
df13= 119.999513288462 df24= 166.666248808928
J: Matrix([[119.999513288462, 0], [0, 166.666248808928]])
J的逆矩阵: Matrix([[0.00833336713288280, 0], [0, 0.00600001504291631]])
F： Matrix([[-0.00140907480732722], [0.00848465542944155]])
Q的变量 Matrix([[1.17423376871538e-5], [-5.09080602106108e-5]])
20 10.0000213478611 9.99997865213886 8.00000032699220 2.00002102086894 11.9999996730078
max△q = 0.0000509080602106108 , max△f = 0.00848465542944155
第12次迭代：
f1= -0.00254539902284989 f2= 0.000587100032618082
df13= 120.000504173862 df24= 166.666927451893
J: Matrix([[120.000504173862, 0], [0, 166.666927451893]])
J的逆矩阵: Matrix([[0.00833329832140666, 0], [0, 0.00599999061174654]])
F： Matrix([[-0.00254539902284989], [0.000587100032618082]])
Q的变量 Matrix([[2.12115694044252e-5], [-3.52259468386458e-6]])
20 10.0000001362917 9.99999986370826 7.99999680439752 2.00000333189422 12.0000031956025
max△q = 0.0000212115694044252 , max△f = 0.00254539902284989
第13次迭代：
f1= -0.000176135271786393 f2= 0.00106058574266399
df13= 120.000083161064 df24= 166.666464434353
J: Matrix([[120.000083161064, 0], [0, 166.666464434353]])
J的逆矩阵: Matrix([[0.00833332755826346, 0], [0, 0.00600000728037212]])
F： Matrix([[-0.000176135271786393], [0.00106058574266399]])
Q的变量 Matrix([[1.46779291435977e-6], [-6.36352217744280e-6]])
20 9.99999866849883 10.0000013315012 7.99999044087534 2.00000822762349 12.0000095591247
max△q = 0.00000636352217744280 , max△f = 0.00106058574266399
第14次迭代：
f1= -0.000318176937362580 f2= 7.33895880102864e-5
df13= 120.000207022088 df24= 166.666549264986
J: Matrix([[120.000207022088, 0], [0, 166.666549264986]])
J的逆矩阵: Matrix([[0.00833331895682423, 0], [0, 0.00600000422646348]])
F： Matrix([[-0.000318176937362580], [7.33895880102864e-5]])
Q的变量 Matrix([[2.65146990374786e-6], [-4.40337838240132e-7]])
20 9.99999866849883 10.0000013315012 7.99999044087534 2.00000822762349 12.0000095591247
max△q = 0.00000265146990374786 , max△f = 0.000318176937362580
1
2
3
4
5
6
7
8
9
10
11
12
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
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126

例3

有回路附加阻力的Barczyk法迭代过程程序

#Barczyk，有回路附加阻力
from sympy import *

q1,r1 = 20,0.375
q2,r2 = 11.5,1.5
q3,r3 = 8.5,2
q4,r4 = 10.5,
1
2
3
4
5
6