《1 前言》
1 前言
湍流模型在求解工程中的复杂水流现象做出了重要贡献, 工程中最常用的就是k -ε 双方程湍流模型, Launder, Rodi, Demuren, Ye Jian等做了很多有意义的工作
[1 ,2 ,3 ,4 ] 。在很多情况下, k -ε 湍流模型能够得出比较符合实际的结果, 但由于使用了Boussinesq的线性本构关系及标量的粘性系数概念, 无法考虑离心力与浮力在不同方向的作用, 因而用于模拟环流、强漩流与浮力分层流等各向异性的湍流时, 误差很大, 甚至定性失真。而在各相异性的湍流模拟中, 代数应力湍流模式精度较高, 而复杂程度又远小于应力通量湍流模型, 所以在近年来得到较快发展, 如Gibson、倪浩清、Ross和Rahman等均给出了较成功的算例
[5 ,6 ,7 ,8 ] , 但这些算例都是针对规则边界的流场及标量场, 对于复杂边界的天然河流及海湾, 则受到很大限制。曲线坐标下环境水动力模型是研究污染物在天然水体扩散输移规律及在水体中的浓度分布的重要工具。
笔者建立了曲线坐标下平面二维水流计算和污染物扩散输移的代数应力湍流模型, 并用于实验室连续弯道污染物扩散输移的数值模拟, 计算了流场及岸边和中心污染物排放的浓度分布, 对该模型的计算值与笔者所建k -ε 模型的计算值以及实测值进行了比较
[9 ] , 结果表明代数应力湍流数学模型明显优于基于各向同性假设条件下的k -ε 模型。
《2 数学模型》
2 数学模型
《2.1通用方程》
2.1通用方程
笛卡儿坐标下深度平均的代数应力湍流模型的通用微分方程为
∂ ∂ t ( H ρ Φ ) + ∂ ∂ x ( H ρ u Φ ) + ∂ ∂ y ( H ρ v Φ ) = ∂ ∂ x ( H Γ Φ ∂ Φ ∂ x ) + ∂ ∂ y ( H Γ Φ ∂ Φ ∂ y ) + S Φ ( 1 ) ∂ ∂ t ( Η ρ Φ ) + ∂ ∂ x ( Η ρ u Φ ) + ∂ ∂ y ( Η ρ v Φ ) = ∂ ∂ x ( Η Γ Φ ∂ Φ ∂ x ) + ∂ ∂ y ( Η Γ Φ ∂ Φ ∂ y ) + S Φ ( 1 )
方程式 (1) 转换到非正交曲线坐标 (ξ , η ) 下
[10 ] , 仅在对流项中使用流速的逆变分量, 而在其他项中使用原始变量, 这样既简化了方程, 又使所有方程仍可写为曲线坐标下的通用方程, 模型的通用微分方程可写为如下通用形式:
∂ ∂ t ( H ρ Φ ) + 1 J ∂ ∂ ξ ( H ρ U Φ ) + 1 J ∂ ∂ η ( H ρ V Φ ) = 1 J ∂ ∂ ξ ( H J Γ Φ α Φ ξ − H J Γ Φ β Φ η ) + 1 J ∂ ∂ η ( H J Γ Φ γ Φ η − H J Γ Φ β Φ ξ ) + S Φ ( ξ ‚ η ) ( 2 ) ∂ ∂ t ( Η ρ Φ ) + 1 J ∂ ∂ ξ ( Η ρ U Φ ) + 1 J ∂ ∂ η ( Η ρ V Φ ) = 1 J ∂ ∂ ξ ( Η J Γ Φ α Φ ξ - Η J Γ Φ β Φ η ) + 1 J ∂ ∂ η ( Η J Γ Φ γ Φ η - Η J Γ Φ β Φ ξ ) + S Φ ( ξ ‚ η ) ( 2 )
式中Φ 为所求问题的因变量;U 和V 分别为Cartesian坐标下x , y 流速u , v 的逆变分量, 仅在对流项中出现;Γ Φα , Γ Φγ , Γ Φβ 为扩散系数;S Φ (ξ , η ) 为源项。模型的控制方程组如表1所示。
表1 通用方程中各变量
Table 1 Variable of general equations
《表1》
方程
Φ
Γ Φα
Γ Φγ
S Φ (ξ , η )
连续
1
0
0
0
x -动量
U
y 2 η Γ ux +x 2 η Γ uy
x 2 ξ Γ uy +y 2 ξ Γ ux
S u
y -动量
v
y 2 η Γ vx +x 2 η Γ vy
x 2 ξ Γ vy +y2 ξ ξ 2 Γ vx
S v
浓度
C
y 2 η Γ Cx +x 2 η Γ Cy
x 2 ξ Γ Cy +y 2 ξ Γ Cx
S C
湍流动能
K
(y 2 η Γ kx +x 2 η Γ ky ) /σ k
(x 2 ξ Γ ky +y 2 ξ Γ kx ) /σ k
S k
耗散率
ε
(y2 η η 2 Γ εx +x2 η η 2 Γ εy ) /σ ε
(x 2 ξ Γ εy +y 2 ξ Γ εx ) /σ ε
S ε
S u = − 1 J [ ρ g H ( y η ∂ z s ∂ ξ − y ξ ∂ z s ∂ η ) ] − τ b x − y η J ∂ ∂ ξ ( ρ H C 1 1 k ) + y ξ J ∂ ∂ η ( ρ H C 1 1 k ) + x ξ J ∂ ∂ η [ H ( μ + μ t y ) J ( y η ∂ v ∂ ξ − y ξ ∂ v ∂ η ) ] − x η J ∂ ∂ ξ [ H ( μ + μ t y ) J ( y η ∂ v ∂ ξ − y ξ ∂ v ∂ η ) ] − 1 J ∂ ∂ ξ { H Γ u β J ∂ u ∂ η } − 1 J ∂ ∂ η { H Γ u β J ∂ u ∂ ξ } ( 3 ) S v = − 1 J [ ρ g H ( x ξ ∂ z s ∂ η − x η ∂ z s ∂ ξ ) ] − τ b y − x ξ J ∂ ∂ η ( ρ H C 1 2 k ) + x η J ∂ ∂ ξ ( ρ H C 1 2 k ) + y η J ∂ ∂ ξ [ H ( μ + μ t x ) J ( x ξ ∂ u ∂ η − x η ∂ u ∂ ξ ) ] − y ξ J ∂ ∂ η [ H ( μ + μ t x ) J ( x ξ ∂ u ∂ η − x η ∂ u ∂ ξ ) ] − 1 J ∂ ∂ ξ ( H Γ v β J ∂ v ∂ η ) − 1 J ∂ ∂ η ( H Γ v β J ∂ v ∂ ξ ) ( 4 ) S k = y η J ∂ ∂ ξ [ H Γ k λ J ( x ξ ∂ k ∂ η − x η ∂ k ∂ ξ ) ] − y ξ J ∂ ∂ η [ H Γ k λ J ( x ξ ∂ k ∂ η − x η ∂ k ∂ ξ ) ] + x ξ J ∂ ∂ η [ H Γ k λ J ( y η ∂ k ∂ ξ − y ξ ∂ k ∂ η ) ] − x η J ∂ ∂ ξ [ H Γ k λ J ( y η ∂ k ∂ ξ − y ξ ∂ k ∂ η ) ] − 1 J ∂ ∂ ξ ( H Γ k β J ∂ k ∂ η ) − 1 J ∂ ∂ η ( H Γ k β J ∂ k ∂ ξ ) + H ρ ( P k − P k v − ε ) ( 5 ) S ε = y η J ∂ ∂ ξ [ H Γ ε λ J ( x ξ ∂ ε ∂ η − x η ∂ ε ∂ ξ ) ] − y ξ J ∂ ∂ η [ H Γ ε λ J ( x ξ ∂ ε ∂ η − x η ∂ ε ∂ ξ ) ] + x ξ J ∂ ∂ η [ H Γ ε λ J ( y η ∂ ε ∂ ξ − y ξ ∂ ε ∂ η ) ] − x η J ∂ ∂ ξ [ H Γ ε λ J ( y η ∂ ε ∂ ξ − y ξ ∂ ε ∂ η ) ] − 1 J ∂ ∂ ξ ( H Γ ε β J ∂ ε ∂ η ) − 1 J ∂ ∂ η ( H Γ ε β J ∂ ε ∂ ξ ) + H ρ [ ε k ( C ε 1 P k − C ε 2 ε ) − P ε v ) ] ( 6 ) S C = S C 0 + 1 J ∂ ∂ ξ { 1 C φ 1 k ε u ′ v ′ ¯ ¯ ¯ ¯ ¯ [ 2 ∂ C ∂ ξ x η y η J − ∂ C ∂ η ( x ξ y η J + x η y ξ J ) ] } + 1 J ∂ ∂ η { 1 C φ 1 k ε u ′ v ′ ¯ ¯ ¯ ¯ ¯ [ 2 ∂ C ∂ η x ξ y ξ J − ∂ C ∂ ξ ( x η y ξ J + x ξ y η J ) ] } + 1 J ∂ ∂ ξ { ( 1 − C φ 2 ) [ u ′ φ ′ ¯ ¯ ¯ ¯ ¯ ¯ ( ∂ u ∂ ξ y 2 η J − ∂ u ∂ η y ξ y η J ) + v ′ φ ′ ¯ ¯ ¯ ¯ ¯ ¯ ( ∂ u ∂ η x ξ y η J − ∂ u ∂ ξ x η y η J ) ] } − 1 J ∂ ∂ η { ( 1 − C φ 2 ) [ u ′ φ ′ ¯ ¯ ¯ ¯ ¯ ¯ ( ∂ u ∂ ξ y ξ y η J − ∂ u ∂ η y 2 ξ J ) + v ′ φ ′ ¯ ¯ ¯ ¯ ¯ ¯ ( ∂ u ∂ η x ξ y ξ J − ∂ u ∂ ξ x η y ξ J ) ] } − 1 J ∂ ∂ ξ { ( 1 − C φ 2 ) [ u ′ φ ′ ¯ ¯ ¯ ¯ ¯ ¯ ( ∂ v ∂ ξ x η y η J − ∂ v ∂ η x η y ξ J ) + v ′ φ ′ ¯ ¯ ¯ ¯ ¯ ¯ ( ∂ v ∂ η x ξ y η J − ∂ v ∂ ξ x 2 η J ) ] } + 1 J ∂ ∂ η { ( 1 − C φ 2 ) [ u ′ φ ′ ¯ ¯ ¯ ¯ ¯ ¯ ( ∂ v ∂ ξ x ξ y η J − ∂ v ∂ η x ξ y ξ J ) + v ′ φ ′ ¯ ¯ ¯ ¯ ¯ ¯ ( ∂ v ∂ η x 2 ξ J − ∂ v ∂ ξ x ξ x η J ) ] } ( 7 ) Γ Φ α = Γ Φ x y 2 η + Γ Φ y x 2 η , Γ Φ γ = Γ Φ y x 2 ξ + Γ Φ x y 2 ξ , Γ Φ β = Γ Φ x y ξ y η + Γ Φ x x ξ x η ( Φ = u , v , C , k , ε ) ‚ Γ u x = 2 μ + 2 μ t 0 + ρ D u , Γ u y = μ + μ t x + ρ D u , Γ v x = μ + 2 μ t y + ρ D v , Γ v y = 2 μ + 2 μ t 0 + ρ D u , Γ C x = μ + ρ u ′ u ′ ¯ ¯ ¯ ¯ ¯ ¯ k / C φ 1 ε , Γ C y = μ + ρ v ′ v ′ ¯ ¯ ¯ ¯ ¯ k / C φ 1 ε , Γ k x = μ + ρ C k u ′ u ′ ¯ ¯ ¯ ¯ ¯ ¯ k / ε , Γ k y = μ + ρ C k v ′ v ′ ¯ ¯ ¯ ¯ ¯ k / ε , Γ ε x = μ + ρ C ε u ′ u ′ ¯ ¯ ¯ ¯ ¯ ¯ k / ε , Γ ε y = μ + ρ C ε v ′ v ′ ¯ ¯ ¯ ¯ ¯ k / ε , Γ k λ = ρ C k u ′ v ′ ¯ ¯ ¯ ¯ ¯ k / ε , Γ ε λ = ρ C ε u ′ v ′ ¯ ¯ ¯ ¯ ¯ k / ε , μ t 0 = ρ C μ 0 k 2 / ε , U = u y η − v x η , V = v x ξ − u x ξ , α = x 2 η + y 2 η , β = x ξ x η + y ξ y η , γ = x 2 ξ + y 2 ξ , J = x ξ y η − y ξ x η , μ t y = ρ C μ 0 k 2 ε ⋅ [ 1 − 2 ( 1 − C 2 ) C 1 k ε ( ∂ u ∂ ξ y η J − ∂ u ∂ η y ξ J ) ] ‚ μ t x = ρ C μ 0 k 2 ε ⋅ [ 1 − 2 ( 1 − C 2 ) C 1 k ε ( ∂ v ∂ η y ξ J − ∂ v ∂ ξ y η J ) ] ‚ P k = − ρ J [ u ′ u ′ ¯ ¯ ¯ ¯ ¯ ¯ ( ∂ u ∂ ξ y η − ∂ u ∂ η y ξ ) + u ′ v ′ ¯ ¯ ¯ ¯ ¯ ( ∂ u ∂ η x ξ − ∂ u ∂ ξ x η ) + u ′ v ′ ¯ ¯ ¯ ¯ ¯ ( ∂ v ∂ ξ y η − ∂ v ∂ η y ξ ) + v ′ v ′ ¯ ¯ ¯ ¯ ¯ ( ∂ v ∂ η x ξ − ∂ v ∂ ξ x η ) ] ‚ P k v = c a ρ U 3 ∗ / H , P ε v = c b ρ U 4 ∗ / H 2 , U ∗ = ( c f ( u 2 + v 2 ) ) 1 / 2 , H = z s − z b , D u = D v = 0 . 2 H U ∗ , c a = c f − 1 / 2 ‚ c b = 3 . 6 C ε 2 c − 0 . 7 5 f C 1 / 2 μ , c f = 0 . 0 0 3 , τ b x = g n 2 u ( u 2 + v 2 ) 1 / 2 / H 1 / 3 , τ b y = g n 2 v ( u 2 + v 2 ) 1 / 2 / H 1 / 3 , S u = - 1 J [ ρ g Η ( y η ∂ z s ∂ ξ - y ξ ∂ z s ∂ η ) ] - τ x b - y η J ∂ ∂ ξ ( ρ Η C 1 1 k ) + y ξ J ∂ ∂ η ( ρ Η C 1 1 k ) + x ξ J ∂ ∂ η [ Η ( μ + μ t y ) J ( y η ∂ v ∂ ξ - y ξ ∂ v ∂ η ) ] - x η J ∂ ∂ ξ [ Η ( μ + μ t y ) J ( y η ∂ v ∂ ξ - y ξ ∂ v ∂ η ) ] - 1 J ∂ ∂ ξ { Η Γ u β J ∂ u ∂ η } - 1 J ∂ ∂ η { Η Γ u β J ∂ u ∂ ξ } ( 3 ) S v = - 1 J [ ρ g Η ( x ξ ∂ z s ∂ η - x η ∂ z s ∂ ξ ) ] - τ y b - x ξ J ∂ ∂ η ( ρ Η C 1 2 k ) + x η J ∂ ∂ ξ ( ρ Η C 1 2 k ) + y η J ∂ ∂ ξ [ Η ( μ + μ t x ) J ( x ξ ∂ u ∂ η - x η ∂ u ∂ ξ ) ] - y ξ J ∂ ∂ η [ Η ( μ + μ t x ) J ( x ξ ∂ u ∂ η - x η ∂ u ∂ ξ ) ] - 1 J ∂ ∂ ξ ( Η Γ v β J ∂ v ∂ η ) - 1 J ∂ ∂ η ( Η Γ v β J ∂ v ∂ ξ ) ( 4 ) S k = y η J ∂ ∂ ξ [ Η Γ k λ J ( x ξ ∂ k ∂ η - x η ∂ k ∂ ξ ) ] - y ξ J ∂ ∂ η [ Η Γ k λ J ( x ξ ∂ k ∂ η - x η ∂ k ∂ ξ ) ] + x ξ J ∂ ∂ η [ Η Γ k λ J ( y η ∂ k ∂ ξ - y ξ ∂ k ∂ η ) ] - x η J ∂ ∂ ξ [ Η Γ k λ J ( y η ∂ k ∂ ξ - y ξ ∂ k ∂ η ) ] - 1 J ∂ ∂ ξ ( Η Γ k β J ∂ k ∂ η ) - 1 J ∂ ∂ η ( Η Γ k β J ∂ k ∂ ξ ) + Η ρ ( Ρ k - Ρ k v - ε ) ( 5 ) S ε = y η J ∂ ∂ ξ [ Η Γ ε λ J ( x ξ ∂ ε ∂ η - x η ∂ ε ∂ ξ ) ] - y ξ J ∂ ∂ η [ Η Γ ε λ J ( x ξ ∂ ε ∂ η - x η ∂ ε ∂ ξ ) ] + x ξ J ∂ ∂ η [ Η Γ ε λ J ( y η ∂ ε ∂ ξ - y ξ ∂ ε ∂ η ) ] - x η J ∂ ∂ ξ [ Η Γ ε λ J ( y η ∂ ε ∂ ξ - y ξ ∂ ε ∂ η ) ] - 1 J ∂ ∂ ξ ( Η Γ ε β J ∂ ε ∂ η ) - 1 J ∂ ∂ η ( Η Γ ε β J ∂ ε ∂ ξ ) + Η ρ [ ε k ( C ε 1 Ρ k - C ε 2 ε ) - Ρ ε v ) ] ( 6 ) S C = S C 0 + 1 J ∂ ∂ ξ { 1 C φ 1 k ε u ′ v ′ ¯ [ 2 ∂ C ∂ ξ x η y η J - ∂ C ∂ η ( x ξ y η J + x η y ξ J ) ] } + 1 J ∂ ∂ η { 1 C φ 1 k ε u ′ v ′ ¯ [ 2 ∂ C ∂ η x ξ y ξ J - ∂ C ∂ ξ ( x η y ξ J + x ξ y η J ) ] } + 1 J ∂ ∂ ξ { ( 1 - C φ 2 ) [ u ′ φ ′ ¯ ( ∂ u ∂ ξ y η 2 J - ∂ u ∂ η y ξ y η J ) + v ′ φ ′ ¯ ( ∂ u ∂ η x ξ y η J - ∂ u ∂ ξ x η y η J ) ] } - 1 J ∂ ∂ η { ( 1 - C φ 2 ) [ u ′ φ ′ ¯ ( ∂ u ∂ ξ y ξ y η J - ∂ u ∂ η y ξ 2 J ) + v ′ φ ′ ¯ ( ∂ u ∂ η x ξ y ξ J - ∂ u ∂ ξ x η y ξ J ) ] } - 1 J ∂ ∂ ξ { ( 1 - C φ 2 ) [ u ′ φ ′ ¯ ( ∂ v ∂ ξ x η y η J - ∂ v ∂ η x η y ξ J ) + v ′ φ ′ ¯ ( ∂ v ∂ η x ξ y η J - ∂ v ∂ ξ x η 2 J ) ] } + 1 J ∂ ∂ η { ( 1 - C φ 2 ) [ u ′ φ ′ ¯ ( ∂ v ∂ ξ x ξ y η J - ∂ v ∂ η x ξ y ξ J ) + v ′ φ ′ ¯ ( ∂ v ∂ η x ξ 2 J - ∂ v ∂ ξ x ξ x η J ) ] } ( 7 ) Γ Φ α = Γ Φ x y η 2 + Γ Φ y x η 2 , Γ Φ γ = Γ Φ y x ξ 2 + Γ Φ x y ξ 2 , Γ Φ β = Γ Φ x y ξ y η + Γ Φ x x ξ x η ( Φ = u , v , C , k , ε ) ‚ Γ u x = 2 μ + 2 μ t 0 + ρ D u , Γ u y = μ + μ t x + ρ D u , Γ v x = μ + 2 μ t y + ρ D v , Γ v y = 2 μ + 2 μ t 0 + ρ D u , Γ C x = μ + ρ u ′ u ′ ¯ k / C φ 1 ε , Γ C y = μ + ρ v ′ v ′ ¯ k / C φ 1 ε , Γ k x = μ + ρ C k u ′ u ′ ¯ k / ε , Γ k y = μ + ρ C k v ′ v ′ ¯ k / ε , Γ ε x = μ + ρ C ε u ′ u ′ ¯ k / ε , Γ ε y = μ + ρ C ε v ′ v ′ ¯ k / ε , Γ k λ = ρ C k u ′ v ′ ¯ k / ε , Γ ε λ = ρ C ε u ′ v ′ ¯ k / ε , μ t 0 = ρ C μ 0 k 2 / ε , U = u y η - v x η , V = v x ξ - u x ξ , α = x η 2 + y η 2 , β = x ξ x η + y ξ y η , γ = x ξ 2 + y ξ 2 , J = x ξ y η - y ξ x η , μ t y = ρ C μ 0 k 2 ε ⋅ [ 1 - 2 ( 1 - C 2 ) C 1 k ε ( ∂ u ∂ ξ y η J - ∂ u ∂ η y ξ J ) ] ‚ μ t x = ρ C μ 0 k 2 ε ⋅ [ 1 - 2 ( 1 - C 2 ) C 1 k ε ( ∂ v ∂ η y ξ J - ∂ v ∂ ξ y η J ) ] ‚ Ρ k = - ρ J [ u ′ u ′ ¯ ( ∂ u ∂ ξ y η - ∂ u ∂ η y ξ ) + u ′ v ′ ¯ ( ∂ u ∂ η x ξ - ∂ u ∂ ξ x η ) + u ′ v ′ ¯ ( ∂ v ∂ ξ y η - ∂ v ∂ η y ξ ) + v ′ v ′ ¯ ( ∂ v ∂ η x ξ - ∂ v ∂ ξ x η ) ] ‚ Ρ k v = c a ρ U * 3 / Η , Ρ ε v = c b ρ U * 4 / Η 2 , U * = ( c f ( u 2 + v 2 ) ) 1 / 2 , Η = z s - z b , D u = D v = 0 . 2 Η U * , c a = c f - 1 / 2 ‚ c b = 3 . 6 C ε 2 c f - 0 . 7 5 C μ 1 / 2 , c f = 0 . 0 0 3 , τ b x = g n 2 u ( u 2 + v 2 ) 1 / 2 / Η 1 / 3 , τ b y = g n 2 v ( u 2 + v 2 ) 1 / 2 / Η 1 / 3 ,
其中, S C0 为污染物源强, k 和ε 分别为深度平均的湍流动能和湍流动能耗散率, H 为水深, z s 为水位, z b 为河床高程, μ 为分子动力粘性系数, − u ′ u ′ ¯ ¯ ¯ ¯ ¯ ¯ ‚ − v ′ v ′ ¯ ¯ ¯ ¯ ¯ , − u ′ v ′ ¯ ¯ ¯ ¯ ¯ - u ′ u ′ ¯ ‚ - v ′ v ′ ¯ , - u ′ v ′ ¯ 分别代表u , v 及交叉方向的Reynolds应力, − u ′ φ ′ ¯ ¯ ¯ ¯ ¯ ¯ ‚ − v ′ φ ′ ¯ ¯ ¯ ¯ ¯ ¯ - u ′ φ ′ ¯ ‚ - v ′ φ ′ ¯ 为方向的通量。
《2.2雷诺应力》
2.2雷诺应力
− u ′ u ′ ¯ ¯ ¯ ¯ ¯ ¯ = 2 C μ 0 k 2 ε ( ∂ u ∂ ξ y η J − ∂ u ∂ η y ξ J ) − C 1 1 k ( 8 ) - u ′ u ′ ¯ = 2 C μ 0 k 2 ε ( ∂ u ∂ ξ y η J - ∂ u ∂ η y ξ J ) - C 1 1 k ( 8 )
− v ′ v ′ ¯ ¯ ¯ ¯ ¯ = 2 C μ 0 k 2 ε ( ∂ v ∂ η x ξ J − ∂ v ∂ ξ x η J ) − C 1 2 k ( 9 ) − u ′ v ′ ¯ ¯ ¯ ¯ ¯ = C μ 0 k 2 ε [ 1 − 2 ( 1 − C 2 ) C 1 k ε ⋅ ( ∂ v ∂ η x ξ J − ∂ v ∂ ξ x η J ) ] ( ∂ u ∂ η x ξ J − ∂ u ∂ ξ x η J ) + C μ 0 k 2 ε [ 1 − 2 ( 1 − C 2 ) C 1 k ε ⋅ ( ∂ u ∂ ξ x η J − ∂ u ∂ η y ξ J ) ] ( ∂ v ∂ ξ y η J − ∂ v ∂ η y ξ J ) ( 1 0 ) u ′ φ ′ ¯ ¯ ¯ ¯ ¯ ¯ = 1 C φ 1 k ε { ( 1 − C φ 2 ) [ u ′ φ ′ ¯ ¯ ¯ ¯ ¯ ¯ ( ∂ u ∂ ξ y η J − ∂ u ∂ η y ξ J ) + v ′ φ ′ ¯ ¯ ¯ ¯ ¯ ¯ ( − ∂ u ∂ ξ x η J + ∂ u ∂ η x ξ J ) ] − u ′ u ′ ¯ ¯ ¯ ¯ ¯ ¯ ( ∂ C ∂ ξ y η J − ∂ C ∂ η y ξ J ) − u ′ v ′ ¯ ¯ ¯ ¯ ¯ ( − ∂ C ∂ ξ x η J + ∂ C ∂ η x ξ J ) } ( 1 1 ) v ′ φ ′ ¯ ¯ ¯ ¯ ¯ ¯ = 1 C φ 1 k ε { ( 1 − C φ 2 ) [ u ′ φ ′ ¯ ¯ ¯ ¯ ¯ ¯ ( ∂ v ∂ ξ y η J − ∂ u ∂ η y ξ J ) + v ′ φ ′ ¯ ¯ ¯ ¯ ¯ ¯ ( − ∂ v ∂ ξ x η J + ∂ u ∂ η x ξ J ) ] − u ′ v ′ ¯ ¯ ¯ ¯ ¯ ( − ∂ C ∂ ξ y η J − ∂ C ∂ η y ξ J ) − v ′ v ′ ¯ ¯ ¯ ¯ ¯ ( − ∂ C ∂ ξ x η J + ∂ C ∂ η x ξ J ) } ( 1 2 ) C μ 0 = 1 3 ( 1 − C 2 ) ( 1 − C 2 + 2 C 1 ) / { C 2 1 + ( 1 − C 2 ) 2 k 2 ε 2 [ ( ∂ u ∂ η x ξ J − ∂ u ∂ ξ x η J ) − ( ∂ v ∂ ξ y η J − ∂ v ∂ η y ξ J ) ] 2 } ( 1 3 ) C 1 1 = 1 3 { ( 1 − C 2 + 2 C 1 ) C 1 + C 1 3 ( ∂ u ∂ η x ξ J − ∂ u ∂ ξ x η J ) ⋅ [ ( ∂ u ∂ η x ξ J − ∂ u ∂ ξ x η J ) − ( ∂ v ∂ ξ y η J − ∂ v ∂ η y ξ J ) ] K ¯ ¯ ¯ 2 ε ¯ 2 } / { C 2 1 + ( 1 − C 2 ) 2 k 2 ε 2 [ ( ∂ u ∂ η x ξ J − ∂ u ∂ ξ x η J ) − ( ∂ v ∂ ξ y η J − ∂ v ∂ η y ξ J ) ] 2 } ( 1 4 ) - v ′ v ′ ¯ = 2 C μ 0 k 2 ε ( ∂ v ∂ η x ξ J - ∂ v ∂ ξ x η J ) - C 1 2 k ( 9 ) - u ′ v ′ ¯ = C μ 0 k 2 ε [ 1 - 2 ( 1 - C 2 ) C 1 k ε ⋅ ( ∂ v ∂ η x ξ J - ∂ v ∂ ξ x η J ) ] ( ∂ u ∂ η x ξ J - ∂ u ∂ ξ x η J ) + C μ 0 k 2 ε [ 1 - 2 ( 1 - C 2 ) C 1 k ε ⋅ ( ∂ u ∂ ξ x η J - ∂ u ∂ η y ξ J ) ] ( ∂ v ∂ ξ y η J - ∂ v ∂ η y ξ J ) ( 1 0 ) u ′ φ ′ ¯ = 1 C φ 1 k ε { ( 1 - C φ 2 ) [ u ′ φ ′ ¯ ( ∂ u ∂ ξ y η J - ∂ u ∂ η y ξ J ) + v ′ φ ′ ¯ ( - ∂ u ∂ ξ x η J + ∂ u ∂ η x ξ J ) ] - u ′ u ′ ¯ ( ∂ C ∂ ξ y η J - ∂ C ∂ η y ξ J ) - u ′ v ′ ¯ ( - ∂ C ∂ ξ x η J + ∂ C ∂ η x ξ J ) } ( 1 1 ) v ′ φ ′ ¯ = 1 C φ 1 k ε { ( 1 - C φ 2 ) [ u ′ φ ′ ¯ ( ∂ v ∂ ξ y η J - ∂ u ∂ η y ξ J ) + v ′ φ ′ ¯ ( - ∂ v ∂ ξ x η J + ∂ u ∂ η x ξ J ) ] - u ′ v ′ ¯ ( - ∂ C ∂ ξ y η J - ∂ C ∂ η y ξ J ) - v ′ v ′ ¯ ( - ∂ C ∂ ξ x η J + ∂ C ∂ η x ξ J ) } ( 1 2 ) C μ 0 = 1 3 ( 1 - C 2 ) ( 1 - C 2 + 2 C 1 ) / { C 1 2 + ( 1 - C 2 ) 2 k 2 ε 2 [ ( ∂ u ∂ η x ξ J - ∂ u ∂ ξ x η J ) - ( ∂ v ∂ ξ y η J - ∂ v ∂ η y ξ J ) ] 2 } ( 1 3 ) C 1 1 = 1 3 { ( 1 - C 2 + 2 C 1 ) C 1 + C 1 3 ( ∂ u ∂ η x ξ J - ∂ u ∂ ξ x η J ) ⋅ [ ( ∂ u ∂ η x ξ J - ∂ u ∂ ξ x η J ) - ( ∂ v ∂ ξ y η J - ∂ v ∂ η y ξ J ) ] Κ ¯ 2 ε ¯ 2 } / { C 1 2 + ( 1 - C 2 ) 2 k 2 ε 2 [ ( ∂ u ∂ η x ξ J - ∂ u ∂ ξ x η J ) - ( ∂ v ∂ ξ y η J - ∂ v ∂ η y ξ J ) ] 2 } ( 1 4 )
C 1 2 = 1 3 { ( 1 − C 2 + 2 C 1 ) C 1 − C 1 3 ( ∂ v ∂ ξ y η J − ∂ v ∂ η y ξ J ) ⋅ [ ( ∂ u ∂ η x ξ J − ∂ u ∂ ξ x η J ) − ( ∂ v ∂ ξ y η J − ∂ v ∂ η y ξ J ) ] k 2 ε 2 } / C 1 2 = 1 3 { ( 1 - C 2 + 2 C 1 ) C 1 - C 1 3 ( ∂ v ∂ ξ y η J - ∂ v ∂ η y ξ J ) ⋅ [ ( ∂ u ∂ η x ξ J - ∂ u ∂ ξ x η J ) - ( ∂ v ∂ ξ y η J - ∂ v ∂ η y ξ J ) ] k 2 ε 2 } /
{ C 2 1 + ( 1 − C 2 ) 2 k 2 ε 2 [ ( ∂ u ∂ η x ξ J − ∂ u ∂ ξ x η J ) − ( ∂ v ∂ ξ y η J − ∂ v ∂ η y ξ J ) ] 2 } ( 1 5 ) C 1 3 = 2 ( 1 − C 2 ) 2 ( 1 − C 2 + 2 C 1 ) / C 1 ( 1 6 ) { C 1 2 + ( 1 - C 2 ) 2 k 2 ε 2 [ ( ∂ u ∂ η x ξ J - ∂ u ∂ ξ x η J ) - ( ∂ v ∂ ξ y η J - ∂ v ∂ η y ξ J ) ] 2 } ( 1 5 ) C 1 3 = 2 ( 1 - C 2 ) 2 ( 1 - C 2 + 2 C 1 ) / C 1 ( 1 6 )
式中 C k =0.24, C 1 =2.2, C 2 =0.55, C ε =0.15, C ε1 =1.44, C ε2 =1.92, C φ1 =3.0, C φ2 =0.5, σ k =1.0, σ k =1.3为湍流经验常数。
《3 模型离散和求解》
3 模型离散和求解
模型的离散采用控制体积法, 控制方程的系数采用乘方格式, 模型求解采用SIMPLEC算法。针对代数应力湍流模型中湍流动能源项容易出现负值而产生发散的问题, 采用文献[11 ] 中的处理办法。
《3.1边界条件》
3.1边界条件
在进口边界, 所有边界条件都按本征条件给出, 即u =u 0 , v =v 0 , k =k 0 , ε =ε 0 。
在出口边界, 给定水位z s , 并设各变量法向导数为零, 即∂ z s ∂ n = ∂ u ∂ n = ∂ v ∂ n = ∂ k ∂ n = ∂ ε ∂ n = 0 ∂ z s ∂ n = ∂ u ∂ n = ∂ v ∂ n = ∂ k ∂ n = ∂ ε ∂ n = 0 。
在壁面边界, 采用壁面函数法
[12 ] , 即用半分析的方法得到的解来近似由壁面到湍流核心区之间的流速、湍流动能、湍流动能耗散率的分布规律, 将壁面的影响 (如壁面应力) 附加到差分方程中 (应先将有关的边界系数置为零) 。
《3.2SIMPLEC算法》
3.2SIMPLEC算法
Step 2 根据进、出口边界的水位, 给定全场水位, 计算全场初始水深H * ;
Step 5 解水深校正方程, 求H ′, 并修正水深H =H * +α h H ′和流速场u , v ;
Step 6 求解各个方向的雷诺应力− u ′ u ′ ¯ ¯ ¯ ¯ ¯ ¯ ‚ − u ′ v ′ ¯ ¯ ¯ ¯ ¯ ‚ − v ′ v ′ ¯ ¯ ¯ ¯ ¯ ‚ − u ′ φ ′ ¯ ¯ ¯ ¯ ¯ ¯ ‚ − v ′ φ ′ ¯ ¯ ¯ ¯ ¯ ¯ - u ′ u ′ ¯ ‚ - u ′ v ′ ¯ ‚ - v ′ v ′ ¯ ‚ - u ′ φ ′ ¯ ‚ - v ′ φ ′ ¯ 及各个方程、各个方向的扩散系数;
Step 7 重复Step 3至Step 6直至收敛;
模型求解中, 为了利于非线性迭代的收敛, 计算中采用亚松弛技术, ADI技术和TDMA算法, 以连续方程的残余质量小于一给定值 (10-5 ) 作为判断收敛的依据。
《4 模型应用》
4 模型应用
模型采用Chang 在实验室连续弯道中进行的一系列水流和污染物扩散输移的试验结果作为验证资料
[13 ] 。弯道断面为矩形并具有光滑底面, 其形状与具体尺寸如图1所示。在试验中, 污染物扩散输移按三种情况进行, 即中心排放, 岸边排放 (左岸排放和右岸排放) , 在第一个弯道进口处排放。水深H =0.115 m, 进口流量Q =0.09486 m3 /s, 污染物初始浓度为零。模型将物理区域划分为166×45个网格, 时间步长取Δt =0.2 s, 糙率取0.015。计算至1 200步, 计算收敛值收敛至规定值。选了4个断面浓度的模型计算值与k -ε 模型计算值以及实测值进行比较
[9 ] , 见图2。
《图1》
图1 实验室连续弯道示意图
Fig.1 Sketch of meandering channel in lab
《图2》
图2 断面浓度比较
Fig.2 Comparison of section concentration
从图2的比较, 可以看出本模型计算的浓度在横向的扩散明显大于基于各向同性假定的k -ε 模型的计算值, 调整Schmidt数σ C 可以加大浓度的横向扩散, 但在没有实测值的情况下, σ C 的调整是没有根据的。连续弯道中污染物浓度的横向扩散还受到弯道二次流的加强, 所以, 平面二维模型在计算曲率比较大的弯曲河流时, 缺陷比较明显。
《5 结论》
5 结论
建立了曲线坐标下平面二维污染物扩散输移的代数应力湍流模型, 以对具有不规则边界水域的湍流各向异性特征进行模拟。将模型用于实验室连续弯道的污染物扩散输移的数值模拟, 分别计算了岸边和中心污染物排放的浓度分布, 对模型的计算值与笔者所建k -ε 湍流模型的计算值以及实测值进行了比较
[9 ] , 结果表明代数应力湍流数学模型明显优于基于各向同性假设条件下的k -ε 模型。由于平面二维模型在计算曲率比较大的弯曲河流时具有比较大的缺陷, 需要发展曲线坐标下三维弯曲河流水流和污染物扩散输移的湍流模型。