有限体积法
基本原理
有限体积法和通量差分在网格划分上基本相同,原理上也差不多。网格划分上,有网格点在网格的中心也有网格点在网格顶点的。不论哪种网格划分,网格边界都是通量点,每一步计算都只有网格点上有值,通量点上的值是利用插值临时计算出来用于推进的。这里我就只做网格点在网格中心的情况。
对于方程
∂
u
∂
t
+
∂
f
∂
x
=
0
\frac{\partial u}{\partial t}+\frac{\partial f}{\partial x}=0
∂t∂u+∂x∂f=0在一个网格单元上做空间积分得
∫
l
r
∂
u
∂
t
+
∂
f
∂
x
d
x
=
0
\int_l^r \frac{\partial u}{\partial t}+\frac{\partial f}{\partial x}dx=0
∫lr∂t∂u+∂x∂fdx=0于是有
∂
∂
t
∫
l
r
u
d
x
+
f
∣
l
r
=
0
\frac{\partial}{\partial t}\int_l^r u dx+f|_l^r=0
∂t∂∫lrudx+f∣lr=0
记
∫
l
r
u
d
x
r
−
l
=
u
ˉ
\frac{\int_l^r udx}{r-l}=\bar u
r−l∫lrudx=uˉ则有
∂
u
ˉ
∂
t
(
r
−
l
)
+
f
∣
l
r
=
0
\frac{\partial \bar u}{\partial t}(r-l)+f|_l^r=0
∂t∂uˉ(r−l)+f∣lr=0
当取
u
ˉ
=
u
i
\bar u=u_i
uˉ=ui且网格单元等长时网格单元的交界面恰好就是
i
−
1
2
i-\frac{1}{2}
i−21和
i
+
1
2
i+\frac{1}{2}
i+21,方程就变成了
∂
u
i
∂
t
+
f
i
+
1
2
−
f
i
−
1
2
Δ
x
=
0
\frac{\partial u_i}{\partial t}+\frac{f_{i+\frac{1}{2}}-f_{i-\frac{1}{2}}}{\Delta x}=0
∂t∂ui+Δxfi+21−fi−21=0
这就是通量差分。
但是有限体积法的思路相比于通量差分更接近真实的物理流动变化,因此这种方法在商业软件中的使用相当广泛。
基本无震荡(ENO)格式
限制器
在前面通量差分格式中可以看到,高阶通量重构的耗散是小于低阶通量耗散的,因此使用高阶通量重构得到的结果在间断位置更加尖锐但是会在间断位置有一个小的尖峰,而使用低阶的通量重构则有较大的耗散,在间断位置就会变得平滑,尖锐位置就被抹平。由于间断处通常是容易发散的位置,所以我们希望间断处的耗散尽可能大一些,而非间断位置耗散尽可能小,不要影响计算结果。这时我们就考虑使用限制器来选择通量的重构方式。根据计算点前后的导数变化情况选择合适的通量构造函数
u
i
n
+
1
=
u
i
n
−
Δ
t
Δ
x
(
φ
(
r
1
)
f
i
+
1
2
(
1
)
−
φ
(
r
2
)
f
i
−
1
2
(
1
)
+
φ
(
r
1
)
f
i
+
1
2
(
2
)
−
φ
(
r
2
)
f
i
−
1
2
(
2
)
)
u_i^{n+1}=u_i^n-\frac{\Delta t}{\Delta x}(\varphi(r_1)f_{i+\frac{1}{2}}^{(1)}-\varphi(r_2)f_{i-\frac{1}{2}}^{(1)}+\varphi(r_1)f_{i+\frac{1}{2}}^{(2)}-\varphi(r_2)f_{i-\frac{1}{2}}^{(2)})
uin+1=uin−ΔxΔt(φ(r1)fi+21(1)−φ(r2)fi−21(1)+φ(r1)fi+21(2)−φ(r2)fi−21(2))
式中
r
r
r表示当地的导数变化情况,通过构造相应的函数
φ
\varphi
φ就可以根据当地导数变化情况选用合适的通量构造方式,从而抑制震荡。
例如当接近间断位置通常函数的导数就会剧烈变化,这时使用高耗散的方法构造通量,而导数变化不剧烈的位置就使用低耗散的方法构造通量,这样就就可以实现抑制间断震荡的同时保证光滑段的计算精度足够高。
基本原理
ENO格式的原理和限制器方法不同,但是思路上却是异曲同工的。首先确定阶数或者说通量点重构需要的求解点数量。以三点格式为例,在构造 f i + 1 2 f_{i+\frac{1}{2}} fi+21时首先必然需要 x i x_i xi,由此得到零阶的近似。然后以 x i x_i xi为中心向两侧加上一个点构造一阶近似,这时有两种构造方法,即使用 x i , x i + 1 x_i,x_{i+1} xi,xi+1或者使用 x i , x i − 1 x_i,x_{i-1} xi,xi−1。这两种方法构造时的会得到两个不同的差商,哪个小就使用哪个,这样可以保证构造的通量始终是导数小的那个,即变化不剧烈的构造方法。再扩展一个点则有更多种构造方法,同样选择差商小的那个来构造最终的通量。通过这种尝试和比较的方法的可以使得构造出来的通量总是最光滑的那一个,最终计算也就不会出现明显的导数剧烈变化,所以这个格式称为基本无震荡(Essentially Non-oscillatory)格式。
由于这两种方法我接触的比较少因此很多具体的东西我也不是很了解,只是了解一个大概,甚至可能这里面有一些错误,所以这里也就不展示解方程的结果了。