學習自F. Moukalled, L. Mangani, M. Darwish所著The Finite Volume Method in Computational Fluid Dynamics - An Advanced Introduction with OpenFOAM and Matlab html
前面幾章講解的關於變量
ϕ
\phi
ϕ web
因爲本章內容繁雜，篇幅較長，故分紅了四部分來說解，各部分主要內容分別爲：交錯網格、同位網格、邊界條件、SIMPLE家族算法。算法
這裏是第一部分，主要講解不可壓縮流動問題的求解在常規網格上所碰到的問題，以及，交錯網格上的SIMPLE算法是如何解決該問題的。網絡
前面幾章所處理的通用守恆方程，能夠改寫成與連續方程或動量方程類似的形式。然而，目前爲止所呈現的數值技術，並不足以求解Navier-Stokes方程，由於求解通常流動的算法要可以處理壓力速度的耦合問題。爲理解該問題，先把連續方程和動量方程寫出來
∂
ρ
∂
t
+
∇
⋅
(
ρ
v
)
=
0
∂
∂
t
(
ρ
v
)
+
∇
⋅
(
ρ
v
v
)
=
−
∇
p
+
∇
⋅
{
μ
[
∇
v
+
(
∇
v
)
T
]
}
+
f
b
\begin{aligned} \frac{\partial\rho}{\partial t}+\nabla\cdot(\rho\bold v)&=0 \\ \frac{\partial}{\partial t}(\rho \bold v)+\nabla\cdot(\rho\bold v\bold v)&= -\nabla p+\nabla\cdot\left\{ \mu\left[ \nabla\bold v+(\nabla\bold v)^\text{T} \right] \right\} + \bold f_b \end{aligned}
∂ t ∂ ρ + ∇ ⋅ ( ρ v ) ∂ t ∂ ( ρ v ) + ∇ ⋅ ( ρ v v ) = 0 = − ∇ p + ∇ ⋅ { μ [ ∇ v + ( ∇ v ) T ] } + f b
(
∇
v
)
T
(\nabla\bold v)^\text{T}
( ∇ v ) T app
從這倆方程中能夠看出，縱然速度場能夠直接使用動量方程來計算，然而出如今動量方程中的壓力場沒法直接從連續方程中計算。該耦合關係（壓力和速度的耦合）是如此強烈卻又如此隱含，爲便於理解，可把方程組寫成以下矩陣形式
A
u
=
[
F
B
T
B
0
]
[
v
p
]
=
[
f
b
0
]
\bold A \bold u=\begin{bmatrix} \bold F & \bold B^\text T \\ \bold B & \bold 0 \end{bmatrix} \begin{bmatrix} \bold v \\ \bold p \end{bmatrix}= \begin{bmatrix} \bold f_b \\ \bold 0 \end{bmatrix}
A u = [ F B B T 0 ] [ v p ] = [ f b 0 ] svg
一種處理方法是，經過把矩陣
A
\bold A
A
L
\bold L
L
U
\bold U
U
A
=
[
F
B
T
B
0
]
=
[
F
0
B
−
B
F
−
1
B
T
]
[
I
F
−
1
B
T
0
I
]
=
L
U
\bold A = \begin{bmatrix} \bold F & \bold B^\text T \\ \bold B & \bold 0 \end{bmatrix}= \begin{bmatrix} \bold F & \bold 0 \\ \bold B & -\bold B\bold F^{-1}\bold B^\text T \end{bmatrix} \begin{bmatrix} \bold I & \bold F^{-1}\bold B^\text T \\ \bold 0 & \bold I \end{bmatrix}= \bold L \bold U
A = [ F B B T 0 ] = [ F B 0 − B F − 1 B T ] [ I 0 F − 1 B T I ] = L U
−
B
F
−
1
B
T
-\bold B\bold F^{-1}\bold B^\text T
− B F − 1 B T 函數
這正是下面要講的迭代求解Navier-Stokes方程的方法的本質，Patankar和Spalding將該技術嵌入在經典的分離式SIMPLE（Semi Implicit Method for Pressure Linked Equations）算法中。學習
求解流程是，將Navier-Stokes方程從新整理成動量方程和壓力方程的形式，而後作離散和順序求解。壓力方程是經過聯合半離散動量方程和連續方程（Schur補矩陣）來建立的。ui
算法是由Picard類型的迭代流程驅動的，其率先求解動量方程，使用的是前次迭代的壓力場。求出的速度場是知足動量方程的可是未必知足質量方程。接下來用該速度場來構建壓力方程，並用所得解來修正壓力和速度場以便知足質量守恆。而後開始新的迭代過程，並重覆上述流程，直到速度和壓力場既知足質量守恆又知足動量守恆爲止。spa
該算法可用以下的矩陣形式來描述
[
I
D
−
1
B
T
0
I
]
[
v
p
]
=
[
v
∗
p
∗
]
\begin{bmatrix} \bold I & \bold D^{-1}\bold B^\text T \\ \bold 0 & \bold I \end{bmatrix}\begin{bmatrix} \bold v \\ \bold p \end{bmatrix}= \begin{bmatrix} \bold v^* \\ \bold p^* \end{bmatrix}
[ I 0 D − 1 B T I ] [ v p ] = [ v ∗ p ∗ ]
[
F
0
B
−
B
D
−
1
B
T
]
[
v
∗
p
∗
]
=
[
f
b
0
]
\begin{bmatrix} \bold F & \bold 0 \\ \bold B & -\bold B\bold D^{-1}\bold B^\text T \end{bmatrix} \begin{bmatrix} \bold v^* \\ \bold p^* \end{bmatrix}= \begin{bmatrix} \bold f_b \\ \bold 0 \end{bmatrix}
[ F B 0 − B D − 1 B T ] [ v ∗ p ∗ ] = [ f b 0 ]
F
−
1
\bold F^{-1}
F − 1
D
−
1
\bold D^{-1}
D − 1
∗
^*
∗
求解：
F
v
∗
=
f
b
\bold F \bold v^*=\bold f_b
F v ∗ = f b
求解：
−
B
D
−
1
B
T
p
∗
=
−
B
v
∗
-\bold B \bold D^{-1} \bold B^\text T \bold p^*=-\bold B \bold v^*
− B D − 1 B T p ∗ = − B v ∗
更新：
v
=
v
∗
−
D
−
1
B
T
p
∗
\bold v=\bold v^*-\bold D^{-1}\bold B^\text T\bold p^*
v = v ∗ − D − 1 B T p ∗
更新：
p
=
p
∗
\bold p=\bold p^*
p = p ∗
這種類型的分裂與SIMPLE系列算法中所用的分裂是相同的，這正是本章的主題。
針對不可壓縮流動問題發展求解算法的困難之處，可經過對上圖所示的一維空間均勻網格問題的離散處理來加深理解。爲簡便起見，假設流動是定常的。簡化的連續和動量方程（寫成守恆形式）爲
∂
(
ρ
u
)
∂
x
=
0
∂
(
ρ
u
u
)
∂
x
=
∂
∂
x
(
μ
∂
u
∂
x
)
−
∂
p
∂
x
\begin{aligned} \frac{\partial (\rho u)}{\partial x}&=0 \\ ~ \\ \frac{\partial(\rho u u)}{\partial x}&=\frac{\partial}{\partial x} \left(\mu \frac{\partial u}{\partial x} \right) -\frac{\partial p}{\partial x} \end{aligned}
∂ x ∂ ( ρ u ) ∂ x ∂ ( ρ u u ) = 0 = ∂ x ∂ ( μ ∂ x ∂ u ) − ∂ x ∂ p
動量方程的離散，首先要對其在單元
C
C
C
∫
V
C
∂
(
ρ
u
u
)
∂
x
d
V
=
∫
V
C
∂
∂
x
(
μ
∂
u
∂
x
)
d
V
−
∫
V
C
∂
p
∂
x
d
V
\int_{V_C}\frac{\partial(\rho u u)}{\partial x}dV= \int_{V_C}\frac{\partial}{\partial x} \left(\mu \frac{\partial u}{\partial x} \right)dV- \int_{V_C}\frac{\partial p}{\partial x}dV
∫ V C ∂ x ∂ ( ρ u u ) d V = ∫ V C ∂ x ∂ ( μ ∂ x ∂ u ) d V − ∫ V C ∂ x ∂ p d V
∫
∂
V
C
(
ρ
u
u
)
d
y
=
∫
∂
V
C
μ
∂
u
∂
x
d
y
−
∫
V
C
∂
p
∂
x
d
V
\int_{\partial V_C}(\rho u u)dy=\int_{\partial V_C}\mu \frac{\partial u}{\partial x}dy - \int_{V_C}\frac{\partial p}{\partial x}dV
∫ ∂ V C ( ρ u u ) d y = ∫ ∂ V C μ ∂ x ∂ u d y − ∫ V C ∂ x ∂ p d V
(
ρ
u
Δ
y
)
e
⏟
m
˙
e
u
e
+
−
(
ρ
u
Δ
y
)
w
⏟
m
˙
w
u
w
=
(
μ
∂
u
∂
x
Δ
y
)
e
−
(
μ
∂
u
∂
x
Δ
y
)
w
−
∫
V
C
∂
p
∂
x
d
V
\underbrace{(\rho u \Delta y)_e}_{\dot m_e} u_e+ \underbrace{-(\rho u \Delta y)_w}_{\dot m_w} u_w= \left( \mu \frac{\partial u}{\partial x}\Delta y \right)_e- \left( \mu \frac{\partial u}{\partial x}\Delta y \right)_w- \int_{V_C}\frac{\partial p}{\partial x}dV
m ˙ e
( ρ u Δ y ) e u e + m ˙ w
− ( ρ u Δ y ) w u w = ( μ ∂ x ∂ u Δ y ) e − ( μ ∂ x ∂ u Δ y ) w − ∫ V C ∂ x ∂ p d V
m
˙
e
u
e
+
m
˙
w
u
w
⏟
C
o
n
v
e
c
t
i
o
n
−
[
(
μ
∂
u
∂
x
Δ
y
)
e
−
(
μ
∂
u
∂
x
Δ
y
)
w
]
⏟
D
i
f
f
u
s
i
o
n
=
−
∫
V
C
∂
p
∂
x
d
V
\underbrace{\dot m_e u_e+\dot m_w u_w}_{Convection}- \underbrace{\left[\left( \mu \frac{\partial u}{\partial x}\Delta y \right)_e- \left( \mu \frac{\partial u}{\partial x}\Delta y \right)_w\right]}_{Diffusion}=- \int_{V_C}\frac{\partial p}{\partial x}dV
C o n v e c t i o n
m ˙ e u e + m ˙ w u w − D i f f u s i o n
[ ( μ ∂ x ∂ u Δ y ) e − ( μ ∂ x ∂ u Δ y ) w ] = − ∫ V C ∂ x ∂ p d V
a
C
u
u
C
+
∑
F
∼
N
B
(
C
)
(
a
F
u
u
F
)
=
b
C
u
−
∫
V
C
∂
p
∂
x
d
V
a_C^u u_C+\sum_{F\sim NB(C)}(a_F^u u_F)=b_C^u-\int_{V_C}\frac{\partial p}{\partial x}dV
a C u u C + F ∼ N B ( C ) ∑ ( a F u u F ) = b C u − ∫ V C ∂ x ∂ p d V
連續方程的離散，也是要對其在單元
C
C
C
∫
V
C
∂
(
ρ
u
)
∂
x
d
V
=
0
\int_{V_C}\frac{\partial(\rho u)}{\partial x}dV=0
∫ V C ∂ x ∂ ( ρ u ) d V = 0
∑
f
∼
n
b
(
C
)
(
ρ
u
Δ
y
)
f
=
(
ρ
u
Δ
y
)
e
−
(
ρ
u
Δ
y
)
w
=
0
\sum_{f\sim nb(C)}(\rho u \Delta y)_f=(\rho u \Delta y)_e - (\rho u \Delta y)_w=0
f ∼ n b ( C ) ∑ ( ρ u Δ y ) f = ( ρ u Δ y ) e − ( ρ u Δ y ) w = 0
∑
f
∼
n
b
(
C
)
m
˙
f
=
m
˙
e
+
m
˙
w
=
0
\sum_{f\sim nb(C)}\dot m_f=\dot m_e + \dot m_w=0
f ∼ n b ( C ) ∑ m ˙ f = m ˙ e + m ˙ w = 0
壓力項離散可由下面兩種方法中的任意一個來實現。第一種方法，經過單點Gauss積分來計算體積分
∫
V
C
∂
p
∂
x
d
V
=
(
∂
p
∂
x
)
C
V
C
\int_{V_C}\frac{\partial p}{\partial x}dV=\left(\frac{\partial p}{\partial x}\right)_CV_C
∫ V C ∂ x ∂ p d V = ( ∂ x ∂ p ) C V C
∫
V
C
∂
p
∂
x
d
V
=
p
E
−
p
W
2
Δ
x
V
C
\int_{V_C}\frac{\partial p}{\partial x}dV=\frac{p_E-p_W}{2\Delta x} V_C
∫ V C ∂ x ∂ p d V = 2 Δ x p E − p W V C
∫
V
C
∂
p
∂
x
d
V
=
∫
∂
V
C
p
d
y
\int_{V_C}\frac{\partial p}{\partial x}dV=\int_{\partial V_C}p dy
∫ V C ∂ x ∂ p d V = ∫ ∂ V C p d y
∫
V
C
∂
p
∂
x
d
V
=
∫
∂
V
C
p
d
y
=
p
e
Δ
y
e
−
p
w
Δ
y
w
=
(
p
e
−
p
w
)
Δ
y
=
(
p
e
−
p
w
)
V
C
Δ
x
\int_{V_C}\frac{\partial p}{\partial x}dV=\int_{\partial V_C}p dy= p_e\Delta y_e - p_w \Delta y_w=(p_e-p_w)\Delta y= (p_e-p_w)\frac{V_C}{\Delta x}
∫ V C ∂ x ∂ p d V = ∫ ∂ V C p d y = p e Δ y e − p w Δ y w = ( p e − p w ) Δ y = ( p e − p w ) Δ x V C
∫
V
C
∂
p
∂
x
d
V
=
[
1
2
(
p
E
+
p
C
)
−
1
2
(
p
C
+
p
W
)
]
V
C
Δ
x
=
p
E
−
p
W
2
Δ
x
V
C
\int_{V_C}\frac{\partial p}{\partial x}dV= \left[\frac12(p_E+p_C)-\frac12(p_C+p_W)\right]\frac{V_C}{\Delta x}= \frac{p_E-p_W}{2\Delta x} V_C
∫ V C ∂ x ∂ p d V = [ 2 1 ( p E + p C ) − 2 1 ( p C + p W ) ] Δ x V C = 2 Δ x p E − p W V C
E
E
E
W
W
W
E
E
E
W
W
W
C
C
C
採用一樣的方法，使用線性插值廓線，並注意到密度是常數，且有
(
Δ
y
)
e
=
(
Δ
y
)
w
=
(
Δ
y
)
C
(\Delta y)_e=(\Delta y)_w=(\Delta y)_C
( Δ y ) e = ( Δ y ) w = ( Δ y ) C
u
E
−
u
W
=
0
u_E-u_W=0
u E − u W = 0
在上上式中，單元
C
C
C
針對上圖所示的壓力和速度值，在點
W
W
W
C
C
C
E
E
E
∫
V
W
∂
p
∂
x
d
V
=
(
p
C
−
p
W
W
)
V
W
2
Δ
x
W
=
(
10
−
10
)
V
W
2
Δ
x
W
=
0
∫
V
C
∂
p
∂
x
d
V
=
(
p
E
−
p
W
)
V
C
2
Δ
x
C
=
(
−
100
+
100
)
V
C
2
Δ
x
C
=
0
∫
V
E
∂
p
∂
x
d
V
=
(
p
E
E
−
p
C
)
V
E
2
Δ
x
E
=
(
10
−
10
)
V
E
2
Δ
x
E
=
0
\begin{aligned} &\int_{V_W}\frac{\partial p}{\partial x}dV=(p_C-p_{WW})\frac{V_W}{2\Delta x_W}=(10-10)\frac{V_W}{2\Delta x_W}=0 \\ &\int_{V_C}\frac{\partial p}{\partial x}dV=(p_E-p_{W})\frac{V_C}{2\Delta x_C}=(-100+100)\frac{V_C}{2\Delta x_C}=0 \\ &\int_{V_E}\frac{\partial p}{\partial x}dV=(p_{EE}-p_{C})\frac{V_E}{2\Delta x_E}=(10-10)\frac{V_E}{2\Delta x_E}=0 \\ \end{aligned}
∫ V W ∂ x ∂ p d V = ( p C − p W W ) 2 Δ x W V W = ( 1 0 − 1 0 ) 2 Δ x W V W = 0 ∫ V C ∂ x ∂ p d V = ( p E − p W ) 2 Δ x C V C = ( − 1 0 0 + 1 0 0 ) 2 Δ x C V C = 0 ∫ V E ∂ x ∂ p d V = ( p E E − p C ) 2 Δ x E V E = ( 1 0 − 1 0 ) 2 Δ x E V E = 0
∫
V
W
∂
u
∂
x
d
V
=
(
u
C
−
u
W
W
)
V
W
2
Δ
x
W
=
(
1
−
1
)
V
W
2
Δ
x
W
=
0
∫
V
C
∂
u
∂
x
d
V
=
(
u
E
−
u
W
)
V
C
2
Δ
x
C
=
(
10
−
10
)
V
C
2
Δ
x
C
=
0
∫
V
E
∂
u
∂
x
d
V
=
(
u
E
E
−
u
C
)
V
E
2
Δ
x
E
=
(
1
−
1
)
V
E
2
Δ
x
E
=
0
\begin{aligned} & \int_{V_W}\frac{\partial u}{\partial x}dV=(u_C-u_{WW})\frac{V_W}{2\Delta x_W}=(1-1)\frac{V_W}{2\Delta x_W}=0 \\ & \int_{V_C}\frac{\partial u}{\partial x}dV=(u_E-u_{W})\frac{V_C}{2\Delta x_C}=(10-10)\frac{V_C}{2\Delta x_C}=0 \\ & \int_{V_E}\frac{\partial u}{\partial x}dV=(u_{EE}-u_{C})\frac{V_E}{2\Delta x_E}=(1-1)\frac{V_E}{2\Delta x_E}=0 \\ \end{aligned}
∫ V W ∂ x ∂ u d V = ( u C − u W W ) 2 Δ x W V W = ( 1 − 1 ) 2 Δ x W V W = 0 ∫ V C ∂ x ∂ u d V = ( u E − u W ) 2 Δ x C V C = ( 1 0 − 1 0 ) 2 Δ x C V C = 0 ∫ V E ∂ x ∂ u d V = ( u E E − u C ) 2 Δ x E V E = ( 1 − 1 ) 2 Δ x E V E = 0
前面公式中所呈現的問題是壓力和速度場的解耦問題。能夠經過把不一樣的變量存儲在交錯位置上來增強它們的耦合特性，由於這樣一來，就不須要用插值手段來計算動量方程中的壓力梯度以及連續方程中的速度場了。這樣的交錯網格以下圖所示，在交錯網格中速度場是存儲在單元面上的（圖a），而壓力場則是存儲在單元形心上的（圖b）。
單元
C
C
C
∑
f
∼
n
b
(
C
)
m
˙
f
=
m
˙
e
+
m
˙
w
=
0
\sum_{f\sim nb(C)}\dot m_f=\dot m_e+\dot m_w=0
f ∼ n b ( C ) ∑ m ˙ f = m ˙ e + m ˙ w = 0
u
e
−
u
w
=
0
u_e-u_w=0
u e − u w = 0
e
e
e
w
w
w
e
e
e
a
e
u
u
e
+
∑
f
∼
N
B
(
e
)
a
f
u
u
f
=
b
e
u
−
V
e
(
∇
p
)
e
=
b
e
u
−
V
e
p
E
−
p
C
δ
x
e
a_e^u u_e+\sum_{f\sim NB(e)}a_f^u u_f=b_e^u-V_e(\nabla p)_e= b_e^u-V_e\frac{p_E-p_C}{\delta x_e}
a e u u e + f ∼ N B ( e ) ∑ a f u u f = b e u − V e ( ∇ p ) e = b e u − V e δ x e p E − p C
接下來要呈現的推導是基於Patankar和Spalding的工做，他們發展了SIMPLE（Semi Implicit Method for Pressure Linked Equations）算法的最初實用形式。
由上圖所示網格上的離散連續方程和動量方程出發
∑
f
∼
n
b
(
C
)
m
˙
f
=
0
a
e
u
u
e
+
∑
f
∼
N
B
(
e
)
a
f
u
u
f
=
b
e
u
−
V
e
(
∂
p
∂
x
)
e
\begin{aligned} \sum_{f\sim nb(C)}\dot m_f &= 0 \\ a_e^u u_e+\sum_{f\sim NB(e)}a_f^u u_f&= b_e^u-V_e\left(\frac{\partial p}{\partial x}\right)_e \end{aligned}
f ∼ n b ( C ) ∑ m ˙ f a e u u e + f ∼ N B ( e ) ∑ a f u u f = 0 = b e u − V e ( ∂ x ∂ p ) e
(
n
)
^{(n)}
( n )
u
(
n
)
u^{(n)}
u ( n )
p
(
n
)
p^{(n)}
p ( n )
∗
^*
∗
a
e
u
u
e
∗
+
∑
f
∼
N
B
(
e
)
a
f
u
u
f
∗
=
b
e
u
−
V
e
(
∂
p
(
n
)
∂
x
)
e
a_e^u u_e^*+\sum_{f\sim NB(e)}a_f^u u_f^*= b_e^u-V_e\left(\frac{\partial p^{(n)}}{\partial x}\right)_e
a e u u e ∗ + f ∼ N B ( e ) ∑ a f u u f ∗ = b e u − V e ( ∂ x ∂ p ( n ) ) e
u
∗
u^*
u ∗
將修正場用上標撇來表示，即，（
u
′
u'
u ′
p
′
p'
p ′
u
=
u
∗
+
u
′
p
=
p
∗
+
p
′
u=u^*+u' \\ p=p^*+p'
u = u ∗ + u ′ p = p ∗ + p ′
m
˙
f
=
m
˙
f
∗
+
ρ
u
′
S
f
x
=
m
˙
f
∗
+
m
f
′
\dot m_f=\dot m_f^*+\rho u' S_f^x=\dot m_f^*+m_f'
m ˙ f = m ˙ f ∗ + ρ u ′ S f x = m ˙ f ∗ + m f ′
m
˙
e
+
m
˙
w
=
m
˙
e
∗
+
m
˙
e
′
+
m
˙
w
∗
+
m
˙
w
′
=
0
\dot m_e+\dot m_w=\dot m_e^*+\dot m_e'+\dot m_w^*+\dot m_w'=0
m ˙ e + m ˙ w = m ˙ e ∗ + m ˙ e ′ + m ˙ w ∗ + m ˙ w ′ = 0
m
˙
e
′
+
m
˙
w
′
=
−
m
˙
e
∗
−
m
˙
w
∗
\dot m_e'+\dot m_w'=-\dot m_e^*-\dot m_w^*
m ˙ e ′ + m ˙ w ′ = − m ˙ e ∗ − m ˙ w ∗
m
˙
e
∗
=
ρ
v
e
∗
⋅
S
e
=
ρ
u
e
∗
S
e
x
=
ρ
u
e
∗
Δ
y
e
m
˙
w
∗
=
ρ
v
w
∗
⋅
S
w
=
ρ
u
w
∗
S
w
x
=
−
ρ
u
w
∗
Δ
y
w
\begin{aligned} \dot m_e^* &= \rho \bold v_e^* \cdot \bold S_e=\rho u_e^* S_e^x = \rho u_e^* \Delta y_e \\ \dot m_w^* &= \rho \bold v_w^* \cdot \bold S_w=\rho u_w^* S_w^x = -\rho u_w^* \Delta y_w \end{aligned}
m ˙ e ∗ m ˙ w ∗ = ρ v e ∗ ⋅ S e = ρ u e ∗ S e x = ρ u e ∗ Δ y e = ρ v w ∗ ⋅ S w = ρ u w ∗ S w x = − ρ u w ∗ Δ y w
m
˙
e
′
=
ρ
v
e
′
⋅
S
e
=
ρ
u
e
′
S
e
x
=
ρ
u
e
′
Δ
y
e
m
˙
w
′
=
ρ
v
w
′
⋅
S
w
=
ρ
u
w
′
S
w
x
=
−
ρ
u
w
′
Δ
y
w
\begin{aligned} \dot m_e' &= \rho \bold v_e' \cdot \bold S_e=\rho u_e' S_e^x = \rho u_e' \Delta y_e \\ \dot m_w' &= \rho \bold v_w' \cdot \bold S_w=\rho u_w' S_w^x = -\rho u_w' \Delta y_w \end{aligned}
m ˙ e ′ m ˙ w ′ = ρ v e ′ ⋅ S e = ρ u e ′ S e x = ρ u e ′ Δ y e = ρ v w ′ ⋅ S w = ρ u w ′ S w x = − ρ u w ′ Δ y w
S
e
x
=
Δ
y
e
S_e^x=\Delta y_e
S e x = Δ y e
S
w
x
=
−
Δ
y
w
S_w^x=-\Delta y_w
S w x = − Δ y w
壓力場並無出如今修正形式的連續方程
m
˙
e
′
+
m
˙
w
′
=
−
m
˙
e
∗
−
m
˙
w
∗
\dot m_e'+\dot m_w'=-\dot m_e^*-\dot m_w^*
m ˙ e ′ + m ˙ w ′ = − m ˙ e ∗ − m ˙ w ∗
a
e
u
u
e
+
∑
f
∼
N
B
(
e
)
a
f
u
u
f
=
b
e
u
−
V
e
(
∂
p
∂
x
)
e
a_e^u u_e+\displaystyle\sum_{f\sim NB(e)}a_f^u u_f=b_e^u-V_e\left(\dfrac{\partial p}{\partial x}\right)_e
a e u u e + f ∼ N B ( e ) ∑ a f u u f = b e u − V e ( ∂ x ∂ p ) e
u
e
+
H
e
(
u
e
)
=
B
e
u
−
D
e
u
(
∂
p
∂
x
)
e
u_e+H_e(u_e) = B_e^u - D_e^u\left(\frac{\partial p}{\partial x}\right)_e
u e + H e ( u e ) = B e u − D e u ( ∂ x ∂ p ) e
H
e
(
u
e
)
=
∑
f
∼
N
B
(
e
)
a
f
u
a
e
u
u
f
B
e
u
=
b
e
u
a
e
u
D
e
u
=
V
e
a
e
u
H_e(u_e) =\sum_{f\sim NB(e)}\frac{a_f^u}{a_e^u} u_f \quad\quad B_e^u = \frac{b_e^u}{a_e^u} \quad\quad D_e^u = \frac{V_e}{a_e^u}
H e ( u e ) = f ∼ N B ( e ) ∑ a e u a f u u f B e u = a e u b e u D e u = a e u V e
u
e
∗
u_e^*
u e ∗
u
e
∗
+
H
e
(
u
e
∗
)
=
B
e
u
−
D
e
u
(
∂
p
(
n
)
∂
x
)
e
u_e^*+H_e(u_e^*) = B_e^u - D_e^u\left(\frac{\partial p^{(n)}}{\partial x}\right)_e
u e ∗ + H e ( u e ∗ ) = B e u − D e u ( ∂ x ∂ p ( n ) ) e
u
e
+
H
e
(
u
e
)
=
B
e
u
−
D
e
u
(
∂
p
∂
x
)
e
u_e+H_e(u_e) = B_e^u - D_e^u\left(\frac{\partial p}{\partial x}\right)_e
u e + H e ( u e ) = B e u − D e u ( ∂ x ∂ p ) e
u
e
∗
+
H
e
(
u
e
∗
)
=
B
e
u
−
D
e
u
(
∂
p
(
n
)
∂
x
)
e
u_e^*+H_e(u_e^*) = B_e^u - D_e^u\left(\frac{\partial p^{(n)}}{\partial x}\right)_e
u e ∗ + H e ( u e ∗ ) = B e u − D e u ( ∂ x ∂ p ( n ) ) e
u
e
′
+
H
e
(
u
′
)
‾
=
−
D
e
u
(
∂
p
′
∂
x
)
e
u_e'+\underline{H_e(u')} = - D_e^u\left(\frac{\partial p'}{\partial x}\right)_e
u e ′ + H e ( u ′ ) = − D e u ( ∂ x ∂ p ′ ) e
w
w
w
u
w
′
+
H
w
(
u
′
)
‾
=
−
D
w
u
(
∂
p
′
∂
x
)
w
u_w'+\underline{H_w(u')} = - D_w^u\left(\frac{\partial p'}{\partial x}\right)_w
u w ′ + H w ( u ′ ) = − D w u ( ∂ x ∂ p ′ ) w
m
˙
e
′
=
ρ
u
e
′
Δ
y
e
\dot m_e' = \rho u_e' \Delta y_e
m ˙ e ′ = ρ u e ′ Δ y e
m
˙
w
′
=
−
ρ
u
w
′
Δ
y
w
\dot m_w' = -\rho u_w' \Delta y_w
m ˙ w ′ = − ρ u w ′ Δ y w
m
˙
e
′
+
m
˙
w
′
=
−
m
˙
e
∗
−
m
˙
w
∗
\dot m_e'+\dot m_w'=-\dot m_e^*-\dot m_w^*
m ˙ e ′ + m ˙ w ′ = − m ˙ e ∗ − m ˙ w ∗
ρ
e
u
e
′
Δ
y
e
+
(
−
ρ
w
u
w
′
Δ
y
w
)
=
−
(
m
˙
e
∗
+
m
˙
w
∗
)
\rho_e u_e' \Delta y_e+(-\rho_w u_w' \Delta y_w)=-(\dot m_e^*+\dot m_w^*)
ρ e u e ′ Δ y e + ( − ρ w u w ′ Δ y w ) = − ( m ˙ e ∗ + m ˙ w ∗ )
u
e
′
u_e'
u e ′
u
w
′
u_w'
u w ′
ρ
e
[
−
H
e
(
u
′
)
‾
−
D
e
u
(
∂
p
′
∂
x
)
e
]
Δ
y
e
−
ρ
w
[
−
H
w
(
u
′
)
‾
−
D
w
u
(
∂
p
′
∂
x
)
w
]
Δ
y
w
=
−
(
m
˙
e
∗
+
m
˙
w
∗
)
\begin{aligned} &\rho_e \left[ -\underline{H_e(u')} - D_e^u\left(\frac{\partial p'}{\partial x}\right)_e \right] \Delta y_e \\ &-\rho_w \left[ -\underline{H_w(u')} - D_w^u\left(\frac{\partial p'}{\partial x}\right)_w \right] \Delta y_w=-(\dot m_e^*+\dot m_w^*) \end{aligned}
ρ e [ − H e ( u ′ ) − D e u ( ∂ x ∂ p ′ ) e ] Δ y e − ρ w [ − H w ( u ′ ) − D w u ( ∂ x ∂ p ′ ) w ] Δ y w = − ( m ˙ e ∗ + m ˙ w ∗ )
ρ
e
[
−
H
e
(
u
′
)
‾
−
D
e
u
(
p
E
′
−
p
C
′
Δ
x
)
e
]
Δ
y
e
+
ρ
w
[
−
H
w
(
u
′
)
‾
−
D
w
u
(
p
C
′
−
p
W
′
Δ
x
)
w
]
(
−
Δ
y
w
)
=
−
(
m
˙
e
∗
+
m
˙
w
∗
)
\begin{aligned} &\rho_e \left[ -\underline{H_e(u')} - D_e^u\left(\frac{p_E'-p_C'}{\Delta x}\right)_e \right] \Delta y_e \\ &+\rho_w \left[ -\underline{H_w(u')} - D_w^u\left(\frac{p_C'-p_W'}{\Delta x}\right)_w \right] (-\Delta y_w)=-(\dot m_e^*+\dot m_w^*) \end{aligned}
ρ e [ − H e ( u ′ ) − D e u ( Δ x p E ′ − p C ′ ) e ] Δ y e + ρ w [ − H w ( u ′ ) − D w u ( Δ x p C ′ − p W ′ ) w ] ( − Δ y w ) = − ( m ˙ e ∗ + m ˙ w ∗ )
−
ρ
e
D
e
u
(
Δ
y
e
Δ
x
e
)
(
p
E
′
−
p
C
′
)
−
ρ
w
D
w
u
(
−
Δ
y
w
Δ
x
w
)
(
p
C
′
−
p
W
′
)
=
−
(
m
˙
e
∗
+
m
˙
w
∗
)
+
ρ
e
H
e
(
u
′
)
‾
Δ
y
e
+
ρ
w
H
w
(
u
′
)
‾
(
−
Δ
y
w
)
\begin{aligned} & -\rho_e D_e^u\left(\frac{\Delta y_e}{\Delta x_e}\right)(p_E'-p_C') -\rho_w D_w^u\left(-\frac{\Delta y_w}{\Delta x_w}\right) (p_C'-p_W')\\ &=-(\dot m_e^*+\dot m_w^*)+\rho_e \underline{H_e(u')}\Delta y_e+\rho_w \underline{H_w(u')}(-\Delta y_w) \end{aligned}
− ρ e D e u ( Δ x e Δ y e ) ( p E ′ − p C ′ ) − ρ w D w u ( − Δ x w Δ y w ) ( p C ′ − p W ′ ) = − ( m ˙ e ∗ + m ˙ w ∗ ) + ρ e H e ( u ′ ) Δ y e + ρ w H w ( u ′ ) ( − Δ y w )
a
C
p
′
p
C
′
+
a
E
p
′
p
E
′
+
a
W
p
′
p
W
′
=
b
C
p
′
a_C^{p'}p_C'+a_E^{p'}p_E'+a_W^{p'}p_W'=b_C^{p'}
a C p ′ p C ′ + a E p ′ p E ′ + a W p ′ p W ′ = b C p ′
a
E
p
′
=
−
ρ
e
D
e
u
Δ
y
e
Δ
x
e
a
W
p
′
=
−
ρ
w
D
w
u
Δ
y
w
Δ
x
w
a
C
p
′
=
−
(
a
E
p
′
+
a
W
p
′
)
b
C
p
′
=
−
(
m
˙
e
∗
+
m
˙
w
∗
)
+
[
ρ
e
Δ
y
e
H
e
(
u
′
)
−
ρ
w
Δ
y
w
H
w
(
u
′
)
]
‾
\begin{aligned} a_E^{p'} &= -\frac{\rho_e D_e^u\Delta y_e}{\Delta x_e} \\ a_W^{p'} &= -\frac{\rho_w D_w^u\Delta y_w}{\Delta x_w}\\ a_C^{p'} &= -(a_E^{p'}+a_W^{p'}) \\ b_C^{p'} &= -(\dot m_e^*+\dot m_w^*)+ \underline{[\rho_e \Delta y_e H_e(u')-\rho_w \Delta y_w H_w(u')]} \end{aligned}
a E p ′ a W p ′ a C p ′ b C p ′ = − Δ x e ρ e D e u Δ y e = − Δ x w ρ w D w u Δ y w = − ( a E p ′ + a W p ′ ) = − ( m ˙ e ∗ + m ˙ w ∗ ) + [ ρ e Δ y e H e ( u ′ ) − ρ w Δ y w H w ( u ′ ) ]
Δ
y
\Delta y
Δ y
1
1
1
使用動量方程和壓力修正方程，能夠得到流動問題的解。在SIMPLE算法中的求解過程是，在每步迭代中，交替地生成壓力和速度場，並前後知足動量方程和連續方程，而後再逐步逼近最終解（連續方程和動量方程這兩個方程都知足的最終解）。該流程並不是是同步進行的，這種求解方程的法子常被稱爲分離式方法（segregated approach）（實際上還有一種把壓力和速度直接耦合起來同步求解的法子呢，不過比較耗內存了）。該分離式的SIMPLE算法的計算流程彙總以下：
從某個猜想的壓力場
p
(
n
)
p^{(n)}
p ( n )
u
(
n
)
u^{(n)}
u ( n )
求解動量方程來得到新的速度場
u
f
∗
u_f^*
u f ∗
a
e
u
u
e
∗
+
∑
f
∼
N
B
(
e
)
a
f
u
u
f
∗
=
b
e
u
−
V
e
(
∂
p
(
n
)
∂
x
)
e
a
w
u
u
w
∗
+
∑
f
∼
N
B
(
e
)
a
f
u
u
f
∗
=
b
w
u
−
V
w
(
∂
p
(
n
)
∂
x
)
w
a_e^u u_e^*+\sum_{f\sim NB(e)}a_f^u u_f^*= b_e^u-V_e\left(\frac{\partial p^{(n)}}{\partial x}\right)_e \\ a_w^u u_w^*+\sum_{f\sim NB(e)}a_f^u u_f^*= b_w^u-V_w\left(\frac{\partial p^{(n)}}{\partial x}\right)_w
a e u u e ∗ + f ∼ N B ( e ) ∑ a f u u f ∗ = b e u − V e ( ∂ x ∂ p ( n ) ) e a w u u w ∗ + f ∼ N B ( e ) ∑ a f u u f ∗ = b w u − V w ( ∂ x ∂ p ( n ) ) w
使用知足動量方程的速度場
u
f
∗
u_f^*
u f ∗
m
˙
f
∗
\dot m_f^*
m ˙ f ∗
m
˙
f
∗
=
ρ
v
f
∗
⋅
S
f
=
ρ
f
u
f
∗
S
f
x
\dot m_f^* = \rho \bold v_f^* \cdot \bold S_f=\rho_f u_f^* S_f^x
m ˙ f ∗ = ρ v f ∗ ⋅ S f = ρ f u f ∗ S f x
使用新的質量流量
m
˙
f
∗
\dot m_f^*
m ˙ f ∗
p
′
p'
p ′
a
C
p
′
p
C
′
+
a
E
p
′
p
E
′
+
a
W
p
′
p
W
′
=
b
C
p
′
a_C^{p'}p_C'+a_E^{p'}p_E'+a_W^{p'}p_W'=b_C^{p'}
a C p ′ p C ′ + a E p ′ p E ′ + a W p ′ p W ′ = b C p ′
更新壓力場和速度場以獲取知足連續方程的場
u
f
∗
∗
=
u
f
∗
+
u
f
′
u
f
′
=
−
D
f
u
(
∂
p
′
∂
x
)
f
p
C
∗
=
p
C
(
n
)
+
p
C
′
m
˙
f
∗
∗
=
m
˙
f
∗
+
m
˙
f
′
m
˙
f
′
=
−
ρ
f
D
f
u
Δ
y
f
(
∂
p
′
∂
x
)
f
\begin{aligned} u_f^{**} &= u_f^* + u_f' \quad\quad u_f'=-D_f^u\left(\frac{\partial p'}{\partial x}\right)_f \\ p_C^* &=p_C^{(n)} + p_C' \\ \dot m_f^{**} &= \dot m_f^* + \dot m_f' \quad\quad \dot m_f'=-\rho_f D_f^u \Delta y_f \left(\frac{\partial p'}{\partial x}\right)_f \end{aligned}
u f ∗ ∗ p C ∗ m ˙ f ∗ ∗ = u f ∗ + u f ′ u f ′ = − D f u ( ∂ x ∂ p ′ ) f = p C ( n ) + p C ′ = m ˙ f ∗ + m ˙ f ′ m ˙ f ′ = − ρ f D f u Δ y f ( ∂ x ∂ p ′ ) f
令
u
(
n
)
=
u
∗
∗
u^{(n)}=u^{**}
u ( n ) = u ∗ ∗
p
(
n
)
=
p
∗
p^{(n)}=p^{*}
p ( n ) = p ∗
返回第2步，並重復該流程直到收斂。
理論是晦澀的，而實例則是生動的，用個簡單的小例子來展現下SIMPLE的強大是再好不過的了。
例1 管道網絡（管網）中的流動
下圖展現的是水力管網系統中的一部分，管道流動中的動量方程能夠寫成
m
˙
=
ρ
u
A
=
−
D
∇
P
\dot m=\rho u A = -D\nabla P
m ˙ = ρ u A = − D ∇ P
D
A
=
0.5
D_A=0.5
D A = 0 . 5
D
B
=
D
F
=
0.4
D_B=D_F=0.4
D B = D F = 0 . 4
D
C
=
D
E
=
0.3
D_C=D_E=0.3
D C = D E = 0 . 3
D
D
=
0.19
D_D=0.19
D D = 0 . 1 9
D
G
=
0.1875
D_G=0.1875
D G = 0 . 1 8 7 5
D
H
=
0.35
D_H=0.35
D H = 0 . 3 5
解
在該系統中，使用質量流量場做爲變量，來替代速度場。這並不會產生什麼問題，由於該例中的動量方程已經忽略掉了對流項和擴散項，僅保留了壓力梯度項，即認爲對流和擴散效應，與壓差驅動力相比很是小，直接忽略了。
使用SIMPLE算法的求解要先用假設的壓力場來求解一個初始速度場，而後再根據這個計算的速度場預測一個壓力場，以知足連續方程。
該流程彙總以下：
由一個猜想的壓力場開始；
使用給定的動量方程計算質量流量；
構造壓力修正方程來知足連續性（質量守恆），並用其修正壓力和速度場。
該例中無需迭代，由於其不含對流項所帶來的非線性效應（若是是用SIMPLE來求解不可壓縮流動，那是要迭代屢次才能收斂的）。
step 1
p
3
(
n
)
=
300
p
6
(
n
)
=
200
p
8
(
n
)
=
120
p_3^{(n)}=300 \quad\quad p_6^{(n)}=200 \quad\quad p_8^{(n)}=120
p 3 ( n ) = 3 0 0 p 6 ( n ) = 2 0 0 p 8 ( n ) = 1 2 0
step 2
m
˙
A
∗
=
D
A
(
p
1
−
p
3
(
n
)
)
=
0.5
∗
(
400
−
300
)
=
50
m
˙
B
∗
=
D
B
(
p
3
(
n
)
−
p
2
)
=
0.4
∗
(
300
−
350
)
=
−
20
m
˙
C
∗
=
D
C
(
p
4
−
p
3
(
n
)
)
=
0.3
∗
(
50
−
300
)
=
−
75
m
˙
D
∗
=
D
D
(
p
3
(
n
)
−
p
6
(
n
)
)
=
0.19
∗
(
300
−
200
)
=
19
m
˙
E
∗
=
D
E
(
p
6
(
n
)
−
p
5
)
=
0.3
∗
(
200
−
300
)
=
−
30
m
˙
F
∗
=
D
F
(
p
7
−
p
6
(
n
)
)
=
0.4
∗
(
80
−
200
)
=
−
48
m
˙
G
∗
=
D
G
(
p
6
(
n
)
−
p
8
(
n
)
)
=
0.1875
∗
(
200
−
120
)
=
15
m
˙
H
∗
=
D
H
(
p
9
−
p
8
(
n
)
)
=
0.35
∗
(
200
−
120
)
=
28
\begin{aligned} & \dot m_A^*=D_A(p_1-p_3^{(n)})=0.5*(400-300)=50 \\ & \dot m_B^*=D_B(p_3^{(n)}-p_2)=0.4*(300-350)=-20 \\ & \dot m_C^*=D_C(p_4-p_3^{(n)})=0.3*(50-300)=-75 \\ & \dot m_D^*=D_D(p_3^{(n)}-p_6^{(n)})=0.19*(300-200)=19 \\ & \dot m_E^*=D_E(p_6^{(n)}-p_5)=0.3*(200-300)=-30 \\ & \dot m_F^*=D_F(p_7-p_6^{(n)})=0.4*(80-200)=-48 \\ & \dot m_G^*=D_G(p_6^{(n)}-p_8^{(n)})=0.1875*(200-120)=15 \\ & \dot m_H^*=D_H(p_9-p_8^{(n)})=0.35*(200-120)=28 \end{aligned}
m ˙ A ∗ = D A ( p 1 − p 3 ( n ) ) = 0 . 5 ∗ ( 4 0 0 − 3 0 0 ) = 5 0 m ˙ B ∗ = D B ( p 3 ( n ) − p 2 ) = 0 . 4 ∗ ( 3 0 0 − 3 5 0 ) = − 2 0 m ˙ C ∗ = D C ( p 4 − p 3 ( n ) ) = 0 . 3 ∗ ( 5 0 − 3 0 0 ) = − 7 5 m ˙ D ∗ = D D ( p 3 ( n ) − p 6 ( n ) ) = 0 . 1 9 ∗ ( 3 0 0 − 2 0 0 ) = 1 9 m ˙ E ∗ = D E ( p 6 ( n ) − p 5 ) = 0 . 3 ∗ ( 2 0 0 − 3 0 0 ) = − 3 0 m ˙ F ∗ = D F ( p 7 − p 6 ( n ) ) = 0 . 4 ∗ ( 8 0 − 2 0 0 ) = − 4 8 m ˙ G ∗ = D G ( p 6 ( n ) − p 8 ( n ) ) = 0 . 1 8 7 5 ∗ ( 2 0 0 − 1 2 0 ) = 1 5 m ˙ H ∗ = D H ( p 9 − p 8 ( n ) ) = 0 . 3 5 ∗ ( 2 0 0 − 1 2 0 ) = 2 8
step 3
∑
∼
k
m
˙
k
∗
\displaystyle\sum_{\sim k}\dot m_k^*
∼ k ∑ m ˙ k ∗
Node 3:
∑
∼
k
(
m
˙
k
∗
)
=
m
˙
B
∗
+
m
˙
D
∗
−
m
˙
A
∗
−
m
˙
C
∗
=
−
20
+
19
−
50
+
75
=
24
Node 6:
∑
∼
k
(
m
˙
k
∗
)
=
m
˙
E
∗
+
m
˙
G
∗
−
m
˙
D
∗
−
m
˙
F
∗
=
−
30
+
15
−
19
+
48
=
14
Node 8:
∑
∼
k
(
m
˙
k
∗
)
=
m
˙
I
∗
−
m
˙
G
∗
−
m
˙
H
∗
=
50
−
15
−
28
=
7
\begin{aligned} & \text{Node 3:\quad}\sum_{\sim k}(\dot m_k^*)=\dot m_B^*+\dot m_D^*-\dot m_A^*-\dot m_C^*=-20+19-50+75=24 \\ & \text{Node 6:\quad}\sum_{\sim k}(\dot m_k^*)=\dot m_E^*+\dot m_G^*-\dot m_D^*-\dot m_F^*=-30+15-19+48=14 \\ & \text{Node 8:\quad}\sum_{\sim k}(\dot m_k^*)=\dot m_I^*-\dot m_G^*-\dot m_H^*=50-15-28=7 \end{aligned}
Node 3: ∼ k ∑ ( m ˙ k ∗ ) = m ˙ B ∗ + m ˙ D ∗ − m ˙ A ∗ − m ˙ C ∗ = − 2 0 + 1 9 − 5 0 + 7 5 = 2 4 Node 6: ∼ k ∑ ( m ˙ k ∗ ) = m ˙ E ∗ + m ˙ G ∗ − m ˙ D ∗ − m ˙ F ∗ = − 3 0 + 1 5 − 1 9 + 4 8 = 1 4 Node 8: ∼ k ∑ ( m ˙ k ∗ ) = m ˙ I ∗ − m ˙ G ∗ − m ˙ H ∗ = 5 0 − 1 5 − 2 8 = 7
∑
∼
k
(
m
˙
k
∗
+
m
˙
k
′
)
=
0
⇒
∑
∼
k
(
m
˙
k
′
)
=
−
∑
∼
k
(
m
˙
k
∗
)
\sum_{\sim k}(\dot m_k^*+\dot m_k')=0 \Rightarrow \sum_{\sim k}(\dot m_k')=-\sum_{\sim k}(\dot m_k^*)
∼ k ∑ ( m ˙ k ∗ + m ˙ k ′ ) = 0 ⇒ ∼ k ∑ ( m ˙ k ′ ) = − ∼ k ∑ ( m ˙ k ∗ )
m
˙
A
′
=
D
A
(
−
p
3
′
)
m
˙
B
′
=
D
B
(
p
3
′
)
m
˙
C
′
=
D
C
(
−
p
3
′
)
m
˙
D
′
=
D
D
(
p
3
′
−
p
6
′
)
m
˙
E
′
=
D
E
(
p
6
′
)
m
˙
F
′
=
D
F
(
−
p
6
′
)
m
˙
G
′
=
D
G
(
p
6
′
−
p
8
′
)
m
˙
H
′
=
D
H
(
−
p
8
′
)
\begin{aligned} & \dot m'_A = D_A(-p'_3) \\ & \dot m_B'=D_B(p_3') \\ & \dot m_C'=D_C(-p_3') \\ & \dot m_D'=D_D(p_3'-p_6') \\ & \dot m_E'=D_E(p_6') \\ & \dot m_F'=D_F(-p_6') \\ & \dot m_G'=D_G(p_6'-p_8') \\ & \dot m_H'=D_H(-p_8') \end{aligned}
m ˙ A ′ = D A ( − p 3 ′ ) m ˙ B ′ = D B ( p 3 ′ ) m ˙ C ′ = D C ( − p 3 ′ ) m ˙ D ′ = D D ( p 3 ′ − p 6 ′ ) m ˙ E ′ = D E ( p 6 ′ ) m ˙ F ′ = D F ( − p 6 ′ ) m ˙ G ′ = D G ( p 6 ′ − p 8 ′ ) m ˙ H ′ = D H ( − p 8 ′ )
p
1
,
2
,
4
,
5
,
7
′
p'_{1,2,4,5,7}
p 1 , 2 , 4 , 5 , 7 ′
節點
3
、
6
、
8
三、六、8
3 、 6 、 8
m
˙
B
′
+
m
˙
D
′
−
m
˙
A
′
−
m
˙
C
′
=
−
(
m
˙
B
∗
+
m
˙
D
∗
−
m
˙
A
∗
−
m
˙
C
∗
)
m
˙
E
′
+
m
˙
G
′
−
m
˙
D
′
−
m
˙
F
′
=
−
(
m
˙
E
∗
+
m
˙
G
∗
−
m
˙
D
∗
−
m
˙
F
∗
)
−
m
˙
G
′
−
m
˙
H
′
=
−
(
m
˙
I
∗
−
m
˙
G
∗
−
m
˙
H
∗
)
\begin{aligned} \dot m_B'+\dot m_D'-\dot m_A'-\dot m_C' &= -(\dot m_B^*+\dot m_D^*-\dot m_A^*-\dot m_C^*) \\ \dot m_E'+\dot m_G'-\dot m_D'-\dot m_F' &= -(\dot m_E^*+\dot m_G^*-\dot m_D^*-\dot m_F^*)\\ -\dot m_G'-\dot m_H' &= -(\dot m_I^*-\dot m_G^*-\dot m_H^*) \end{aligned}
m ˙ B ′ + m ˙ D ′ − m ˙ A ′ − m ˙ C ′ m ˙ E ′ + m ˙ G ′ − m ˙ D ′ − m ˙ F ′ − m ˙ G ′ − m ˙ H ′ = − ( m ˙ B ∗ + m ˙ D ∗ − m ˙ A ∗ − m ˙ C ∗ ) = − ( m ˙ E ∗ + m ˙ G ∗ − m ˙ D ∗ − m ˙ F ∗ ) = − ( m ˙ I ∗ − m ˙ G ∗ − m ˙ H ∗ )
D
B
(
p
3
′
)
+
D
D
(
p
3
′
−
p
6
′
)
−
D
A
(
−
p
3
′
)
−
D
C
(
−
p
3
′
)
=
−
(
m
˙
B
∗
+
m
˙
D
∗
−
m
˙
A
∗
−
m
˙
C
∗
)
D
E
(
p
6
′
)
+
D
G
(
p
6
′
−
p
8
′
)
−
D
D
(
