扩散问题的有限体积法

CV
A
∫ ∫ ∫ ∫ ∂
∂t
(ρφ )dV
CV
+
nr
A
→
• (ρφ U )dA
=
nr
A
•
(Γ
gradφ )dA
+
Sφ dV
CV
扩散问题的有限体积法
◆有限体积法输运方程的物理意义
∫ ∫ ∫ ∫ ∂
∂t
(ρφ )dV
CV
+
nr
A
•
→
(ρφ U )dA
=
nr
A
•
(Γ
gradφ )dA
+
Sφ dV
CV
总的变 化率
+
⎜⎜⎝⎛
Γw Aw
δxWP
⎟⎟⎠⎞φW
+
⎜⎜⎝⎛
Γs As
δy SP
⎟⎟⎠⎞φS
+
⎜⎜⎝⎛
Γn
δy
An
PN
⎟⎟⎠⎞φn
+ Su
aPφP = aEφE + aW φW + aSφS + aNφN + Su
aE
=
Γe Ae δxPE
aW
=
Γ w
Aw
δxWP
aS
=
Γs As δySP
aP = aE + aW + aS + aN − S P
A
CV
一维稳态扩散问题的有限体积法
◆一维稳态纯扩散方程
d ⎜⎛ Γ dφ ⎟⎞ + S = 0
dx ⎝ dx ⎠
◆节点划分（P点）
有限体积法的第一步是把求解域划为离散的控制容积。
一维稳态扩散问题的有限体积法
▼控制容积的取法
方法A：一种是把控制容积的界面放在相邻2个节点中间（先划分节点）
方法B：一种是把控制容积的中心节点放在控制容积的几何中心（先划 分控制容积）
一维非稳态扩散问题的有限体积法
◆一维非稳态扩散问题的控制微分方程
ρc
∂T ∂t
=
∂ ∂x
⎜⎛ ⎝
λ
∂T ∂x
⎟⎞ + S ⎠
◆节点划分
一维非稳态扩散问题的有限体积法
◆一维非稳态问题的控制微分方程
∫ ∫ ∫ ∫ ∫ ∫ t+Δt
t CV
ρc
∂T ∂t
dVdt
=
t +Δt t
CV
∂ ∂x
⎜⎛ ⎝
λ
∂T ∂x
aN
=
Γn An δy PN
三维稳态扩散问题的有限体积法
◆三维稳态纯扩散方程
∂ ∂x
⎜⎛ Γ ⎝
∂φ ∂x
⎟⎞ ⎠
+
∂ ∂y
⎜⎜⎝⎛ Γ
∂φ ∂y
⎟⎟⎠⎞
+
∂ ∂z
⎜⎛ Γ ⎝
∂φ ∂z
⎟⎞ ⎠
+
S
=
0
◆节点划分
三维稳态扩散问题的有限体积法
◆三维稳态纯扩散离散方程
a
φ
P
P
=
a
φ
E
E
+
aW
φ W
+ aSφS
=
Γe Ae
φE −φP δxPE
y方向n，s两个界面
⎜⎜⎝⎛
ΓA
∂φ ∂y
⎟⎟⎠⎞ n
=
Γn An
φ N
−φP
δy NP
⎜⎛ ⎝
ΓA
∂φ ∂x
⎟⎞ ⎠w
=
Γw Aw
φ P
−
φ W
δxWP
⎜⎜⎝⎛
ΓA
∂φ
∂y
⎟⎟⎠⎞ s
=
Γs As
φP −φS δy PS
二维稳态扩散问题的有限体积法
Γe
Ae
φ −φ
b = SΔx
aP
=
a
0 p
= θ (aW
+ aE )
一维非稳态扩散问题的有限体积法
◆显式格式
θ =0
b
=
SΔx
=
Su
+
S
PT
0 p
[ ( )] aPTP
=
a E TE0
+ aW TW0
+
a
0 p
−
aE + aW − S P
TP0 + Su
aP
=
a
0 P
a
0 P
=
ρc
Δx Δt
aW
=
λw δxWP
aW
=
λ w
δxWP
在计算中心节点温度 TP 时， 用到了上一时刻的 TW ，TE ，TP 的值，也 同时用到了当前时刻的 TW ，TE 的值（未知）。因此，它不能直接计 算出结果，必须在每个时刻联立求解所有节点的离散方程才能得到
结果，所以它属于隐式格式。此格式被称为Crank-Nicolson格式， 它是一种半隐格式。
首先将微分方程在控制容积上进行积分，利用高斯定理把体 积分转化为控制容积边界界面上的面积分，然后通过对界面 上的参数的近似而得到最终的离散方程。
对界面上的有关参数的近似方法是确定最终离散格式的核心
◆方程的求解（举例）
在每个节点都建立上述离散（对于内部节点，并不需要在每个 节点上重复上述过程，内部节点的离散方程适用于所有内部节 点，而对边界节点则须重新按上述过程进行推导，因为不同的 边界节点界面上有关参数的近似处理方法不同），得到一个线
◆Crank-Nicolson格式（半隐式格式 ）
θ = 0.5
a P TP
= aE
TE
+ TE0 2
+ aW
TW
+ TW0 2
+
⎡ ⎢⎣a
0 p
− ⎜⎛ aE ⎝
+ aW 2
⎟⎠⎞⎥⎦⎤TP0
+b
aP
=
aW
+ aE 2
+
a
0 P
− SP 2
aE
=
λ e
δxPE
a
0 P
=
ρc
Δx Δt
b
=
Su
+
S PTP0 2
−
(1
−
θ
)
λ e
δxPE
−
(1
−
θ
)
λ w
δxWP
⎥⎦⎤TP0
+ SΔx
[ ] [ ] aPTP = aE θTE + (1−θ )TE0 + aW θTW + (1−θ )TW0 +
[ ] a
0 p
− (1−θ
)aE
− (1−θ
)aW
TP0
+b
aW
=
λw δxWP
a
0 P
=
ρc
Δx Δt
aE
=
λe δxPE
当 λ 为常数，且采用均匀网格时 δxPE = δxWP = Δx
ρc
Δx Δt
−
λ Δx
−
λ Δx
>
0
Δt
<
ρc
Δx 2
2λ
显式格式稳定性条件
当采用显式格式计算时，如果希望采用较小的空间步长以取得更为 精确的结果，则时间步长将非常小。这将使得计算时间很长。因此， 一般不推荐显式格式。
一维非稳态扩散问题的有限体积法
aE
=
λe δxPE
在计算中心节点温度 TP 时，只用到了上一时刻的TW ，TE ，TP 的值， 因此它叫显式格式，可直接由初始温度分布计算出其它时刻的温度 分布 。
一维非稳态扩散问题的有限体积法
◆显式格式稳定性条件
a
0 p
− (aE
+ aW
− SP )>
0
SP < 0
a
0 p
−
(aE
+
aW
)
>
0
一维非稳态扩散问题的有限体积法
◆全隐式格式
θ =1
aPTP
=
aETE
+ aW TW
+
a
T0 0
pP
+ Su
aP
= aW
+ aE
+
a
0 P
− SP
a
0 P
=
ρc
Δx Δt
aW
=
λw δxWP
aE
=
λe δxPE
在计算中心节点温度 TP 时，用到了当前时刻的 TW ，TE 的值（未知）。 因此，它是全隐格式。在每个时刻，必须对所有节点的离散方程同时 求解，才能得到各节点的温度值，给定一个初始值，就可以逐时计算。 该式中所有节点温度的系数都是正值，因此它是无条件稳定的。但它 的精度是一阶（对时间项来说），所以要想提高计算精度，必须采用 较小的步长。全隐式格式一般被推荐作为非稳态问题的格式。
∫ ∫ ∫ ∫ ∫ ∫ ∫ ∫ Δt
∂ ∂t
(ρφ )dVdt
CV
+
Δt
nr
A
•
→
(ρφ U )dAdt
=
Δt
A
nr
•
(Γ
gradφ )dAdt
+
Sφ dVdt
Δt CV
去
掉
对
流
∫ ∫ ∫ ∫ ∫ ∫ Δt
∂ ∂t
(ρφ )dV dt
CV
=
Δt
A
nr
项
• (Γ
gradφ )dAdt
+
Sφ dVdt
Δt CV
+ aNφN
+
a
φ
B
B
+
aT
φ T
+ Su
aE
=
Γe Ae
δxPE
aN
=
Γn An
δy PN
aW
=
Γw Aw
δxWP
aB
=
Γb Ab
δz BP
aS
=
Γs As
δySP
aT
=
Γt At
δz PT
aP = aE + aW + aS + aN + aB + aT − S P
非稳态扩散问题的有限体积法
◆非稳态流动与传热的输运方程最通用的形式的积分方程
一维稳态扩散问题的有限体积法
◆方程的离散
∫ ∫ d ⎜⎛Γ dφ ⎟⎞dV + SdV = ⎜⎛ΓA dφ ⎟⎞ − ⎜⎛ΓA dφ ⎟⎞ + SΔV = 0
CV dx ⎝ dx ⎠
CV
⎝ dx ⎠e ⎝ dx ⎠w
⎜⎛ΓA dφ
⎝ dx
⎟⎞ ⎠e
=
Γ e
Ae
φ E
−φP
δxPE
⎜⎛ ΓA ⎝
dφ
t
权系数 0 ≤ θ ≤ 1
∫ [ ] ( ) IT
=
t +Δt
t TP dt
=
θTP
+
1 − θ TP0
( ) ρc
TP − TP0 Δt
Δx
=
θ
⎜⎜⎝⎛
λe
A
TE − TP
δxPE
−
λ
w
A
TP − TW
δxWP
⎟⎟⎠⎞ +
(1
−
θ
)⎜⎜⎝⎛
λe
A
TEwenku.baidu.com − TP0
δxPE
−
λw
A
TP0 − TW0
⎜⎛ Γ ⎝
∂φ
∂x
⎟⎞ ⎠
+
∂ ∂y
⎜⎜⎝⎛ Γ
∂φ
∂y
⎟⎟⎠⎞
+
S
=
0
◆高斯定理把体积分转换为面积分得
⎜⎛ ⎝
ΓA
∂φ ∂x
⎟⎞ ⎠e
−
⎜⎛ ⎝
ΓA
∂φ ∂x
⎟⎞ ⎠w
+
⎜⎜⎝⎛
ΓA
∂φ ∂y
⎟⎟⎠⎞ n
−
⎜⎜⎝⎛
ΓA
∂φ ∂y
⎟⎟⎠⎞ s
+
SΔV
=
0
x方向e，w两个界面
⎜⎛ ⎝
ΓA
∂φ
∂x
⎟⎞ ⎠e
+
Γw Aw δxWP
− SP
⎟⎟⎠⎞φP
=
⎜⎜⎝⎛
Γe Ae δxPE
⎟⎟⎠⎞φE
+
⎜⎜⎝⎛
Γw Aw δxWP
⎟⎟⎠⎞φW
+ Su
a
φ
P
P
=
a
φ
E
E
+
aW
φ W
+ Su
aP = aE + aW − SP
aW
=
Γw Aw
δxWP
aE
=
Γe Ae
δxPE
一维稳态扩散问题的有限体积法
◆方程离散的步骤
流体力学数值方法
第六讲
扩散问题的有限体积法
扩散问题的有限体积法
◆通用形式流动与传热问题守恒形式的输运方程
∂(ρφ )
∂t
+
→
div(ρ U φ)
=
div(Γ
gradφ )
+
Sφ
瞬变项
对流项
扩散项
源项
▼在应用有限体积法（控制容积法）进行数值求 解时，通常首先将通用公式在一个容积上进行积 分，将微分方程转化为积分方程，然后采用不同 的近似方式在控制容积的边界上对积分项进行处 理，从而得到不同的差分格式。
E
P
δxPE
−
Γw
Aw
φP − φW δxWP
+
Γn
An
φN −φP δy NP
−
Γs
As
φP −φS δy PS
+ SΔV
=0
SΔV
=
Su
+
S
φ
P
P
⎜⎜⎝⎛
Γe Ae
δxPE
+
Γw Aw
δxWP
+
Γs As
δy SP
+
Γn An
δy PN
− S P ⎟⎟⎠⎞φP
=
⎜⎜⎝⎛
Γe Ae
δxPE
⎟⎟⎠⎞φE
性方程组。求解该方程组即可得每个节点的 φ 值。
二维稳态扩散问题的有限体积法
◆二维稳态纯扩散方程
∂ ∂x
⎜⎛ Γ ⎝
∂φ
∂x
⎟⎞ ⎠
+
∂ ∂y
⎜⎜⎝⎛ Γ
∂φ
∂y
⎟⎟⎠⎞
+
S
=
0
◆节点划分
有限体积法的第一步是把求解域划为离散的控制容积。
二维稳态扩散问题的有限体积法
◆控制方程在控制容积上积分
∂ ∂x
一维非稳态扩散问题的有限体积法
◆半隐式格式稳定性条件
为保证计算结果物理上的真实性和有界性，式中各节点温度的系数 须为正
a
0 p
− ⎜⎛ aE ⎝
+ aW 2
⎟⎞ ⎠
>
0
Δt
<
ρc
Δx 2
λ
半隐式格式稳定性条件
Crank-Nicolson格式的稳定性条件与显式格式比，并没有很大的改 善，但此格式采用是中心差分（对时间项），其截差为二阶，它的 精度比显式格式好。
∂ ∂t
(ρφ )dVdt
CV
+
Δt
A
nr
•
→
(ρφ U )dAdt
=
Δt
A
nr
•
(Γ
gradφ )dAdt
+
Sφ dVdt
Δt CV
扩散问题的有限体积法
◆稳态纯扩散
( ) div Γ gradφ + Sφ = 0
∫ div(Γ gradφ)dV + ∫ Sφ dV = 0
CV
CV
∫ nr • (Γ gradφ)dA + ∫ Sφ dV = 0
外法线方 向的对流 通量
物理意义
因对流而 引起的净 减少量
内法线方 向的扩散 通量
物理意义
扩散而引 起的净增 加量
源项引 起的的 增加率
扩散问题的有限体积法
◆稳态输运方程
∫ nr
→
• (ρφ U )dA
=
∫ nr
• (Γ
gradφ)dA +
∫ Sφ dV
A
A
CV
◆非稳态输运方程
∫ ∫ ∫ ∫ ∫ ∫ ∫ ∫ Δt
扩散问题的有限体积法
◆有限体积法求解过程
∫ ∫ ∫ ∫ CV
∂(ρφ )
∂t dV
+
→
div(ρ U φ)dV
CV
=
div(Γ
CV
gradφ )dV
+
Sφ dV
CV
高斯定理
∫ div(ar)dV = ∫ nr • ardA
CV
A
∫
div(
ρ
→
U
φ
)dV
=
∫ nr
→
• (ρφ U )dA
CV
A
∫ div ( Γ grad φ ) dV = ∫ nr • ( Γ grad φ ) dA
dx
⎟⎞ ⎠w
=
Γw Aw
φP − φW δxWP
Γw
=
ΓW
+ ΓP 2
Γe
=
Γ E
+Γ P
2
中心差分格式
S Δ V = S u + S Pφ P
一维稳态扩散问题的有限体积法
◆方程的离散
Γe
Ae
φ E
− φP
δxPE
−
Γw
Aw
φ −φ PW δxWP
+ Su
+ SPφP
=0
⎜⎜⎝⎛
Γe Ae δxPE
δxWP
⎟⎞ ⎟⎠
+
SΔx
一维非稳态扩散问题的有限体积法
◆一维非稳态问题的控制微分方程
[ ] ⎢⎡ρc
⎣
Δx Δt
+
θ
⎜⎜⎝⎛
λ e
δxPE
+
λ w
δxWP
⎟⎟⎠⎞⎥⎦⎤TP
=
λ e
δxPE
θTE
+ (1 − θ )TE0
+
[ ] λw
δxWP
θTW
+ (1 − θ )TW0
+
⎢⎣⎡ρc
Δx Δt
∫ ∫ ( ) CV
⎡ ⎢⎣
t +Δt t
ρc
∂T ∂t
dt ⎥⎦⎤dV
=
ρc
TP
− TP0
ΔV
一维非稳态扩散问题的有限体积法
◆一维非稳态问题的控制微分方程
( ) ∫ ∫ ρc TP −TP0 ΔV
=
t t
+
Δt
⎜⎜⎝⎛
λ e
A
TE − TP
δxPE
−
λ w
A
TP − TW
δxWP
⎟⎟⎠⎞dt
+
t+Δt SΔV dt
⎟⎞dVdt ⎠
+
t +Δt
SdVdt
t CV
∫ ∫ ∫ ∫ ⎡
CV ⎢⎣
t t
+ Δt
ρc
∂T ∂t
⎤ dt ⎥⎦dV
=
t +Δt t
⎢⎣⎡⎜⎝⎛ λA
∂T ∂x
⎟⎞ ⎠e
−
⎜⎛ λA ⎝
∂T ∂x
⎟⎞ ⎠w
⎤ ⎥⎦dt
+
t+Δt SΔVdt
t
∂T
= TP
−
T
0 p
∂t
Δt
Tp0
t时刻的温度
TP
当前 t + Δt 时刻的节点温度