流体力学有限元分析中的边界条件处理

流体力学有限元分析中的边界条件处理

梁启国 尹敏镐

赵永凯 高殿荣

燕山大学

大庆石化总厂

摘 要 阐述了流体力学有限元分析中应用流函数 - 涡量法时典型边界条件的处理方 法 ,并给出计算实例 .

关键词 有限元 涡 - 流函数法 边界条件 分类号 O357 . 1

引言

求解不可压缩粘性流体二维流动问题的数值方法有速度 - 压力法 、涡 - 流函数法和流函数法 ,其中涡 - 流函数法应用较广1 ～4 .

在应用涡 - 流函数法进行有限元分析时 ,数值边界条件的处理不但影响解的精度而且影响解的稳定性2 . 本文结合工程实际给出了在有限元分析中确定几种典型边界条件的方法.

0 边界条件的处理

如图 1 所示的沿背部台阶的不可压缩粘性流 体二维流动. C 点为角点.

在涡 - 流函数方程的求解过程中 , 不仅需要 知道边界上的流函数和涡量值 ,有时还要对边界 上的涡量值不断地进行修正. 下面分别讨论各种 边界上的流函数值和涡量值. 设Δ x 和Δy 分别为

x 方向和 y 方向的有限元划分网格间距 , i 和 j 分别

为 x 方向和 y 方向的步长指标.

1) 壁面边界处理

若壁面是不可渗透的 ,则沿壁面的流函数 Ψ

为常数. 一阶精度的壁涡公式为 :

1 图 1 沿背部台阶的边界条件

2 (Ψ i , j

- Ψ i , j +1 ) Ωi , j

c c

+ O (Δy ) ( 1)

=

(Δy ) 2

c 二阶精度的壁涡公式为 :

3 (Ψ i , j - Ψ i , j +1 ) Ωi , j + 1 + O ( (Δy ) 2 ) Ωi , j c c

c ( 2) = (Δy ) 2

2 c 图 2

用内点上的 Ψ i , j +1 , Ωi , j +1 表示边界上的涡量值 . 二阶壁涡公 c c

式在大网格雷诺数和变网格情况下往往引起数值不稳定和计算不收敛 ; 一阶壁涡公式虽

然精度低 ,但计算稳定 、易于收敛.

另外一种二阶壁涡公式5 为 :

8Ψ i , j +1 - 7Ψ i , j - Ψ i , j +2 + O ( (Δy ) 2 )

Ωi , j c c c

( 3)

=

2 (Δy )

2

c

6 ,

7 M ΩQ + ( K + A )Ω = s

( 4) ( 5)

式中 ,ΩQ = 9

Ω/ 9t , K , M , A , q 和 s 为系数矩阵. 在求解式 (4) 和 (5) 时 ,为确定壁面涡量 , 首先由式 (4) 和 t = t 0 时刻的涡量初值来计算 t = t 0 +Δt 时刻的流函数 ,有 Ψ t 0 +Δt = K - 1 ( M Ωt 0 - q )

( 6) 然后再用式 (4) 和新的流函数值 Ψ t 0 +Δt 求得壁涡

Ωt 0 +Δt = M - 1 ( K Ψ t 0 +Δt + q )

( 7)

最后 ,从式 (5) 出发 ,并利用新的壁涡值可求 t 0 +Δt 时刻其内部各点的涡量值 . 以后每增

加一个时间步长都按上述方法确定壁涡 ,所以也称壁面涡的时间迭代法.

边界 B 6 的处理和 B 1 相同. 2) 对称边界处理

由于流动的对称性 ,图 1 中的边界 B 5 为对称边界 ;而对称线必为流线 ,因此在对称边 界上的流函数 ΨB 为常数. 5

对称边界是滑移边界 ,显然

9 9u

v

v B = 0 = 0

= 0

9x

9y

5

B

B

5

5

由此可求对称边界上的涡量

9 9u

v ΩB = = 0

- 9x 9y

5 B 5

3) 上边界处理

上边界 B 3 一般为以下两种情形 :其一是非滑移固体壁面 ,比如管道流动的上边界就

属于这种情形 . 可采用对 B 1 的处理方法来处理 . 其二是滑移面 ,如飞机在飞行中的无界流体的绕流问题 :此时 y 方向的流场应是无限的 ,但数值计算时却要在有限区域内进行 ,于是在离开绕流物体足够远处确定一个上边界 B 3 , 其上的边界条件应与无穷远条件相匹配 . 滑移面 B 3 可看作滑移固体壁面 ,通常可采用下述两种处理方法 :

第 3 期

赵永凯 等 流体力学有限元分析中的边界条件处理 191

(i ) 在 B 3 上规定为无穷远平行来流 , u = V ∞ , v = 0 , 流函数满足 Ψ i , J 值由流量确定 .

为求 B 3 上的涡量值 ,可推出类似的壁涡公式 : 一阶精度公式为

= co n st . 常数 2 (Ψ i , J -

1 - Ψ i , J + V ∞Δy ) Ωi , J + O (Δy )

= (Δy ) 2 二阶精度公式为

3 (Ψ i , J -

1 - ΨO ( (Δy )

2 ) Ωi , J =

(

. J = V ∞. 因为

= 0 , u i , J = 0 . 可见该条件比 V i , J = 0 条件要宽些. 但对 2

Ψ =

Ωi , Ω

Ψ i , J = Ψ i , J +1 + V ∞Δy

如果流量已知 ,也可以用 Ψ i , J = co n st .

4) 进流边界处理

进流边界 B 2 也称上游边界 . 其边界条件应按进流物理条件确定. 有以下几种情况 : (i ) 规定 u 0 , j = V ∞ , v 0 , j = 0

9u

9v u 0 , j

V ∞ , Ω0 , j = 0 . 因为Ωi , j =(ii ) 规定 , 相当于规定

= -

=

9y

9x

0 , j

0 , j

9v = 0 , 这比规定 v 0 = 0 要放松些.

, j 9x

0 , j

92Ψ

(iii ) 规定 Ψ0 , j , 9

9u v

= 0 . 据此Ω0 =也就给定了. 显

-

= -

, j

9y 2

9x

9y

0 , j

0 , j

0 , j

然 ,此时已不是均匀来流. 通常 Ψ 值由边界层计算而得.

92Ψ

9u

9u

9v

( iv ) 规 定 Ψ0 , j . 而 Ω0 , j 其 中 ,

= - =

,

=

9y 2

9y

9y 9x

0 , j

0 , j

0 , j

0 , j

92Ψ

92Ψ

9v

9v . , = 0 , 可知

= = -

.

9x 2

9x 2

9x 9x

0 , j

0 , j

1 , j

0 , j

1 , j

最终得到

92Ψ

92Ψ

Ω0 , j =

-

9x 2

9y 2

1 , j

0 , j

Ψ

0 , j - 2Ψ1 , j + Ψ2 , j Ψ0 , j +1 - 2Ψ0 , j + Ψ0 , j - 1 = -

-

(Δ x ) 2 (Δy ) 2