基于通量差分裂的圣维南方程组离散方法

向小华,吴晓玲,王春艳,王船海

(1.河海大学水文水资源学院,江苏南京 210098;2.河海大学水文水资源与水利工程科学国家重点实验室,江苏 南京 210098;3.北京中水新华国际工程咨询有限公司,北京 100044)

20世纪70年代以来,计算气体动力学在间断模拟方面取得了突破性进展.20世纪80年代后,Glaister等试图将气体力学中成熟的格式移植到水流计算中[1-3],也取得了一些进展,但所采用的数值方法大多侧重于理论研究,所研究的河道大都相当简单,所论述的方法不能直接应用于实际河道的计算.直到近年来才有一些将现代高精度方法应用于实际河道一维模拟的研究成果[4-6]出现.本文以通量差分裂方法为基础,提出了一种圣维南方程组的高精度离散方法.

1 基本方程

矩形断面棱柱河道圣维南方程组的矢量形式为

其中

式中:Ut——U 对时间t的导数;F(U)x——F(U)对空间 x的导数;h——水深;q——单宽流量;g——重力加速度;n——河段糙率;R——断面水力半径,可以近似用水深 h代替;zb——河底高程.方程(1)为散度形式,即数学上的守恒形式方程.该方程也可以写成非守恒形式

其中

式中:Ux——U对空间x的导数;U——对流速度.

2 数值方法

约定

式中fi,n表示控制体i中tn时刻变量f的值.

首先以简单的线性方程来说明通量差分裂的基本原理.一种波动从tn时刻传播至tn+1时刻的情况如图1所示,采用有限体积法,假设在tn时刻控制体内变量的分布如图1(a)所示,考虑常系数对流方程

式中 ¯u为平均对流速度,常数.

图1 波动传播情况Fig.1 Sketch of wave p ropagation

如果¯u＞0,波动正方向传播,tn→tn+1时刻的传播情况如图1(b)所示,tn+1时刻的波形如图1(d)所示,tn+1时刻控制体i内变量U的值可以采用分段平均的方法进行计算:

其中 :δ Ui-1/2,n=Ui,n-Ui-1,n,可以理解为变量 U 的变差;δ Fi-1/2=¯δ Ui-1/2,n=Fi,n-Fi-1,n.

式(5)可以理解为:在Δt时段内,控制体内平均量等于初始时刻的平均量与波动带来的变差之差,其中变差对平均量的影响与波动速度和时间步长相关.

将式(5)和式(6)相结合则可以写成通用形式

由式(7)的定义可得

式(8)的含义是,可根据波动的方向将界面通量差分裂为2部分,只有波动进入控制体的通量差对控制体具有贡献.从式(7)和式(8)可知,随着时间的推移,控制体内变量平均值的变化完全由时段内进入控制体的通量差确定.

以上述离散方法的原理为基础,离散式(2)所示圣维南方程组,得[7]

式中:Λ——雅可比矩阵A的特征值构成的对角矩阵,其中系数考虑了熵修正;L——高精度限制器.熵修正和限制器的具体表达式见文献[7].

3 算 例

3.1 理想溃坝模型

为了检验本文所述格式的有效性,以常用理想溃坝模型为算例进行计算与分析对比.考察上游水深1.0m、下游水深0.01m时的溃坝情况.多种格式所得溃坝7 s后的水位对比[8-9]如图2所示.总体上看各种格式计算稳定,没有产生振荡现象.从图2可以看出,FOU格式精度最差,间断分辨率最低,另外3种格式都具有较高的精度.计算结果表明:FD-WAF格式精度略低,而MUSCL格式与本文FDS-LW格式精度相当,都很好地捕捉到了间断;本文格式(FDS-LW)的精度更高,间断处接近精确解,优于其他3种格式.

图2 多种格式模拟水位对比Fig.2 Comparison of water levels in different schemes

3.2 水槽试验

Tseng[10]提出了一种地形非常不规则的水槽试验算例,算例所用水槽长1600m,断面为矩形,宽1m,槽底地形很不规则且水槽坡度较大,水槽底部地形(表1)类似于实际河道底部地形,水槽糙率为0.033.

表1 水槽底部高程分布Table 1 Elevation distribution at bottom of flume m

本文模拟所采用的空间步长为7.5m,courant数CFL=0.9.设计的水流条件为恒定流.计算条件如下:初始条件

恒定边界条件

上述条件包括跨临界流等流态.传统的离散格式处理跨临界流存在困难[11-12].图3表明,本文方法所论述的格式在处理跨临界流时可以得到与实测值吻合程度非常好的计算结果.这说明本文方法在处理具有复杂地形的河道水流时也具有较好的适用性.

图3 不规则地形水槽恒定流模拟Fig.3 Irregular channel flow simulated by steady flow in flume

4 结 语

为了提高水流模型的计算精度,本文将气体动力学方法引入水流计算中,并通过离散得出了圣维南方程组的一般离散模型.算例结果表明,本文方法较目前其他一些方法有一定的优势,可以在实际河道水流计算中加以应用.

[1]GLAISTER P.Approximate Riemann solutions of the shallow water equations[J].Journal of Hydraulic Research,1988,26(3):293-306.

[2]TOMAS C R,DELGADO A D,ENRIQUE D F N.A family of stable numerical solvers for the shallow water equations with source terms[J].Comput Methods Appl Mech Engrg,2003,192:203-225.

[3]王志力,耿艳芬,金生.带源项浅水方程的通量平衡离散[J].水科学进展,2005,16(3):373-379.(WANG Zhi-li,GENG Yanfen,JIN Sheng.Flux balance method for shallowwater equationwith source terms[J].Advances in Water Science,2005,16(3):373-379.(in Chinese))

[4]DELIS A I,SKEELS C P.TVD schemes for open channel flow[J].Int J Numer Meth Fluids,1998,26:791-809.

[5]SENKA V,LUKA S.Upwind schemes with exact conservation property for one-dimensional open channel flow equations[J].Society for Industrial and Applied Mathematics,2003,24(5):1630-1649.

[6]BOJAN C,LADO K,NELIDA C.Balanced semi-implicit Roe scheme for open channelflow equations[J].Ann Univ Ferrara,2008,54:245-258.

[7]向小华.基于隐式TVD格式的一、二维耦合流域河网模型研究[D].南京:河海大学,2009.

[8]TORO E F.Shock-cap turing methods free-surface shallow flows[M].London:John Wiley and Sons Ltd.,2001.

[9]ROE P L.Approximate Riemann solvers,parameter vectors,and difference schemes[J].J Comput Phys,1981,43(2):357-372.

[10]TSENG M H.Improved treatment of source terms in TVD scheme for shallow water equations[J].Advances in Water Resources,2004,27:617-629.

[11]FREITAG M,MORTON K W.The Preissmann box scheme and itsmodification for transcritical flows[J].Int JNumer Meth Engng,2007,70:791-811.

[12]茅泽育,相鹏,黄江川.Preissmann四点隐格式对计算混合流动的适定性分析[J].科学技术与工程,2007,7(3):343-347.(MAO Ze-yu,XIANG Peng,HUANG Jiang-chuan.Study on the invalidity of Preissmann scheme for transcritical flow[J].Science Technology and Engineering,2007,7(3):343-347.(in Chinese)).