《1 引言》

1 引言

由于经济的迅速发展, 大量有害于人类和其他生物生息的工农业废水和生活污水未经充分处理而直接排入河流、湖泊、水库、海洋和土壤, 使得地表水和地下水等天然水体受到严重污染。天然水体的河流、湖泊、海洋等都具有一定的自净能力, 如何利用其自净能力, 节约环境保护费用, 是一项非常有意义的研究课题。而弄清污染物在天然水体中的扩散输移规律及其在水体中的浓度分布是研究这一问题的前提, 环境水动力模型是研究这一问题的有效工具。同时, 对弯曲及不规则岸边界条件下的河流、海湾的水动力特性及物质输运规律的研究是环境水力学的基本课题之一。天然河流、海湾的边界曲折、地形复杂, 采用坐标变换是解决问题的途径之一。目前, 在多数N-S方程的坐标变换中, 流速全部采用逆变分量, 这样就增加了方程的复杂程度。于是忽略方程中的非正交项, 利用正交变换下的方程进行数值求解[1,2]。自Patankar和Spalding[3]发展了SIMPLE算法以来, 该方法被广泛应用于不可压缩流场的数值模拟, 而且该方法还得到了进一步的发展, 主要有SIMPLER算法[4]、SIMPLEC算法[5]、SIMPLEX算法[6]和SIMPLET算法[7]等。这些模型均成功地应用于速度—压力耦合的流场计算, 深度平均的浅水流动模型是在静压假定下导出的, 一般流体模型中的速度—压力耦合也就转换成浅水流动模型中的速度—水深耦合[8]

笔者采用非正交曲线坐标下污染物扩散输移的数学模型, 引入污染物的对流扩散方程, 形成非正交曲线坐标系下平面二维水流计算的SIMPLEC算法, 用于实验室连续弯道的水流及污染物扩散输移的数值模拟, 分别计算了流场及岸边和中心污染物排放的浓度场。

《2 数学模型》

2 数学模型

在笛卡儿坐标下深度平均的水流和污染物扩散输移的通用微分方程为

t(ΗρΦ)+x(ΗρuΦ)+y(ΗρvΦ)=x(ΗΓΦΦx)+y(ΗΓΦΦy)+SΦ(1)

方程式式 (1) 通过泊松方程进行坐标变换[2], 转换到曲线坐标 (ξ, η) 下, 仅在对流项中使用流速的逆变分量, 而在其他项中使用原始变量, 这样既简化了方程, 又使所有方程仍可写为非正交曲线坐标下的通用方程, 模型的微分方程可写成如下通用形式:

t(ΗρΦ)+1Jξ(ΗρUΦ)+1Jη(ΗρVΦ)=1Jξ(ΗΓΦJ(αΦξ-βΦη))+1Jη(ΗΓΦJ(-βΦξ+γΦη))+SΦ(ξη)(2)

式中:Φ为所求问题的因变量;UV分别为直角坐标下流速uv的逆变分量, 仅在对流项中出现;ΓΦ为扩散系数;SΦ为源项。当Φ表示某一特定量时, 式 (2) 亦赋予特定的意义, 扩散系数ΓΦ和源项SΦ对此特定量有特定的意义和表达式, 如表1所示。

表1 各控制方程的扩散系数与源项 Table 1 Diffusive coefficient and source term of each control equation

《表1》

方程ΦΓΦSΦ
连续100
动量ξuμeffSu
动量ηνμeffSv
浓度方程CμeffσcSc
湍流动能kμeffσkSk
湍流动能耗散率εμeffσεSε

在式 (2) 及上表中:

Su=-1J[ρgΗ(yηzsξ-yξzsη)]+1Jξ{ΗμeffJ[yη2uξ-yξyηuη]}+1Jξ{ΗμeffJ[yξxηvη-xηyηvξ]}+1Jη{ΗμeffJ[xξyηvξ-xξyξvη]}+1Jη{ΗμeffJ[yξ2uη-yξyηuξ]}-τbξSv=-1J[ρgΗ(xξzsη-xηzsξ)]+1Jξ{ΗμeffJ[xη2vξ-xξxηvη]}+1Jξ{ΗμeffJ[xξyηuη-xηyηuξ]}+1Jη{ΗμeffJ[yξxηuξ-xξyξuη]}+1Jη{ΗμeffJ[xξ2vη-xξxηvξ]}-τbηSk=Η(Gk0+Gkv-ρε)Sε=Η[εk(Cε1Gk0-Cε2ρε)+Gεv]Gk0=2μeff(uξyηJ-uηyξJ)2+2μeff(-vξxηJ+vηxξJ)2+μeff[(vξyηJ-vηyξJ)+(-uξxηJ+uηxξJ)]2,U=uyη-vxηV=vxξ-uyξ,α=xη2+yη2,β=xξxη+yηyξ,γ=xξ2+yξ2,J=xξyη-xηyξ,Gkv=ckρU*3Η,Gεv=cερU*4Η2,U*=cf(u2+v2),Η=zs-zb,ck=1cf,cε=3.6Cε2cf0.75Cμμeff=ρ(v+vt),vt=Cμk2ε,τbξ=gn2Uu2+v2Η1/3τbη=gn2Vu2+v2Η1/3,cf=0.003

上列方程组中, Sc为污染源强, 根据具体的污染物排放形式而定;uv分别为笛卡儿坐标下xy方向的深度平均的流速分量;kε分别为深度平均的湍流动能和湍流动能耗散率;Η为水深;zs为水位;zb为河床高程;v为分子运动粘性系数;vt为漩涡运动粘性系数;μeff为有效粘性系数;Cμ, Cε1 , Cε1, σk, σε, σc为经验常数, 按表2取值。

表2 经验常数 Table 2 Empirical constant

《表2》

CμCε1Cε2σkσεσc
0.091.441.921.01.31.0

《3 模型离散》

3 模型离散

用控制容积法导出离散方程, 控制容积界面上函数的取值及其导数的构造方法决定了对流—扩散方程所采用的格式, 这里采用了混合格式。

aΡΦΡ=aEΦE+aWΦW+aΝΦΝ+aSΦS+a0(3)

其中:

aE=max(|12Fe|,De)-12Fe,aW=max(|12Fw|,Dw)-12Fw,aΝ=max(|12Fn|,Dn)-12Fn,aS=max(|12Fs|,Ds)-12Fs,Fe=(ρUΔη)e,Fw=(ρUΔη)w,Fn=(ρUΔξ)n,Fs=(ρUΔξ)s,De=(αJΓΦΔηΔξ)eDw=(αJΓΦΔηΔξ)wDn=(γJΓΦΔξΔη)nDs=(γJΓΦΔξΔη)sa0=ScJΔξΔη-[(ΓΦJβΦηΔη)we+(ΓΦJβΦξΔξ)sn]+ΔξΔηΔtΗ¯ΡΦ0aΡ=aE+aW+aΝ+aS-SΡ+ΔξΔηΔtΗ¯Ρ

《4 边界条件》

4 边界条件

在进口边界, 所有边界条件都按本征条件给出, 即:

u=u0,v=0,k=k0,ε=ε0,C=C0

在出口边界, 给定出口处的水位zs, 令v=0, 其余变量都取法向导数为零, 即:

zn=un=Cn=kn=εn=0

在壁面边界, k-ε双方程湍流模型是高Reynolds数模型, 只适用于充分发展的湍流, 而在近壁处, 粘性效应起主导作用。另外, 由于近壁处物理量变化很快, 需要非常密的网格, 才能真实地模拟实际流动。因此采用壁面函数法[9], 即用半分析的方法得到的解来近似由壁面到湍流核心区之间的流速、湍动动能、湍动动能耗散率的分布规律, 将壁面的影响 (如壁面应力) 附加到差分方程中 (应先将有关的边界系数置为零) [10]

《5 SIMPLEC算法求解模型》

5 SIMPLEC算法求解模型

SIMPLEC算法的基本步骤如下:

Step 1 计算坐标变换的相关系数;

Step 2 根据进口、出口边界的水位, 给定初始条件及全场水位, 并计算全场初始水深H*;

Step 3 解动量方程求un, vn;

Step 4 解k-ε湍流模型求k, ε;

Step 5 解水深校正方程求H′, 并修正水深H=H*+αhH′;

Step 6 修正流速u, v, 及有效粘性系数μeff, 然后重复step 2 … step 5步直至收敛。

Step 7 求解浓度方程, 直至收敛。

模型求解中, 为了利于非线性迭代的收敛, 计算中采用欠松弛技术。计算中采用ADI技术和TDMA算法, 以连续方程的误差小于给定的值作为判断收敛的依据。在恒定流计算中, 以相邻时间步之间全场各点流速误差的最大绝对值小于给定的值作为达到定常状态的依据。

《6 模型验证》

6 模型验证

模型采用Chang (1971) 在实验室连续弯道中进行的一系列水流和污染物扩散输移的试验结果作为验证资料。弯道断面为矩形并具有光滑底面, 其形状与具体尺寸如图1所示。在其试验中, 污染物扩散输移进行了两种情况, 岸边排放和中心排放, 排放口在第一个弯道进口处, 如图1所示。水深为H=0.115 m, 进口流速U0=0.366 m/s, 污染物初始浓度为零。模型将物理区域划分为167×23个网格, 时间步长取Δt=5 s, 糙率取0.015。计算至800步, 计算收敛值收敛至规定值, 其计算流速场如图2。选了4个断面的流速的实测值与计算值进行比较, 如图3。河流中心污染物排放的断面浓度比较图见图4, 岸边排放的断面浓度断面图见图5。

《图1》

图1 实验室连续弯道示意图 (单位:m) 
Fig.1 Sketch of meandering channel
 in laboratory (Unit: m)

图1 实验室连续弯道示意图 (单位:m) Fig.1 Sketch of meandering channel in laboratory (Unit: m)   

《图2》

图2 实验室连续弯道计算流场
Fig.2 Computed velocity field of meandering
 channel in laboratory

图2 实验室连续弯道计算流场 Fig.2 Computed velocity field of meandering channel in laboratory  

《图3》

图3 断面流速比较
Fig.3 Comparison of section velocity

图3 断面流速比较 Fig.3 Comparison of section velocity   

《图4》

图4 中心排放的断面浓度比较
Fig.4 Comparison of section concentration
 discharged in center

图4 中心排放的断面浓度比较 Fig.4 Comparison of section concentration discharged in center   

《图5》

图5 岸边排放的断面浓度比较
Fig.5 Comparison of section concentration
 discharged from bank

图5 岸边排放的断面浓度比较 Fig.5 Comparison of section concentration discharged from bank  

《7 结论》

7 结论

建立了非正交曲线坐标下水流和污染物扩散输移的数学模型, 对实验室连续弯道中流场、污染物岸边与河中心排放的浓度场进行了数值模拟, 并对计算值与断面实测值进行比较, 结果吻合良好。