[推荐]：二维有限元水流数值模拟

wangwenlong

二维浅水流动的有限元进行数值模拟

江春波，安晓谧，张庆海
(清华大学 水利水电工程系)

摘 要：建立了并行计算系统工作站机群，实现了浅水流动的并行数值模拟。采用信息传递接口(MPI)为通讯库，实现了基于图论的区域剖分方法，并且提出了一种显隐交错的数值格式，该格式具有数值计算精度高，适合并行计算的特点。计算了较低雷诺数下（R<300）方柱绕流的情况。并与Suzuki及Breuer的结果进行了比较。

关键词：并行计算;浅水方程;网格剖分；显隐交错格式

基金项目：国家自然科学基金(59979013)；清华大学基础研究基金(JC2001006)资助项目。
作者简介：江春波(1960-)，男，吉林人，副教授，研究方向：计算水力学。

计算流体力学是应用并行计算技术最活跃的领域之一[1]。在我国，并行计算技术在流体力学领域，对于河道，海洋方面的应用研究还很少。可以说，水力学及河流动力学的并行计算研究在国内基本上还是一个空白。本文首次实现了浅水流动的并行数值模拟，并显示了很好的可扩展性和并行性，为应用这一技术进行大规模洪水预报和水环境的数值模拟进行了探索。

1 控制方程

流动控制方程采用二维浅水流动的连续方程与动量方程：


(1)


(2)


式中：η为水位；H为水深；ui为i轴方向垂向平均流速；g为重力加速度；ν为涡粘性系数；ρ为水密度；τsi为沿i轴方向表面切应力；τbi为沿i轴方向底部切应力，i=1，2.

2 并行计算环境

2.1 硬件环境 采用网络并行机群系统COW（Cluster Of Workstation）。COW以普通PC机作为工作站节点，网络互联形式选择丰富，具有投资风险小、编程方便、系统结构灵活、性能价格比高和可扩充性好等优点。本文的硬件环境为4台工作站和1台100兆的交换机。

2.2 软件环境 并行计算操作系统选择Linux，节点机之间通信采用MPI1.2[2]的国际标准。MPI（Message Passing Interface）是1994年5月发布的一种消息传递接口，它实际上是一个消息传递函数库的标准说明。MPI有如下优点：可移植性和易用性；完备的异步通信功能；有正式和详细的精确定义。并行编程过程还采用了BLAS[3]和LAPACK[4]两个标准的线性代数基本运算函数库。BLAS（Basic Linear Algebra Subprograms）能够在每个计算环境中实现最优化，LAPACK 是一个专门用于高性能计算机的线性代数函数库。应用这些标准的函数库进行模块化的编程将极大的提高计算效率和编程效率，从而使得并行程序的性能大大提高。

3 方柱绕流数值模拟

方柱绕流流动分离点的位置是固定的，不像圆柱绕流流动那样随雷诺数而变化。本文针对Re≤300的方柱绕流进行多种工况模拟，并用Breuer[5]和Suzuki[7]的结果进行验证。

3.1 计算域及边界条件参数 如图1所示，坐标原点取在左下角点处。D为方柱边长；H为平板间距；L为方柱中心和入口处的距离；S为平板总长度。

  
图1 方柱绕流计算域示意 图2 方柱绕流计算网格生成示意

为了使边界不影响方柱附近的流动，取：S=40D， H=8D，L=10D.进口处给定抛物线型流速分布，出口处给定水深为常数。在方柱表面和平板壁面给定固壁边界条件u=v=0。

网格生成采用 Delaunay[7]算法。图2为网格生成的示意图。生成3种网格见表1，综合考虑计算时间和计算精度的要求，实际计算时采用网格SC2.

表1 方柱绕流的3种计算网格


--------------------------------------------------------------------------------

网格编号 节点数 单元数 圆柱附近网格尺度 最大网格尺度

--------------------------------------------------------------------------------

SC1 5640 11008 0.05D 0.4D
SC2 22288 44032 0.01D 0.1D
SC3 88796 173280 0.002D 0.05D

--------------------------------------------------------------------------------


区域分解采用一种基于图论的算法，它是由G.L.Miller[8]等 人在1994年提出的。特点是利用了节点的坐标，将区域分解的最优化问题转化为解线性方程组。这种算法主要基于Radon点的概念，通过点集的投影变换和对中心点的求解，最后实现了剖分。

3.2 计算结果 按流速场流态分以下几种情况验证。

3.2.1 层流阶段 在Re=1～50的情况下，流动属于层流，且为定常，当Re=1,5时，由于流体粘性很大，流速很慢，流经方柱时无分离点产生。随着Re的增加，圆柱后面出现一对稳定的涡，如图3(Re=1,30)所示。

验证时针对回流区长度Lr进行了验证，和Breuer的经验公式Lr/D=-0.065+0.0554Re (5<Re<60) 完全吻合，如图4所示。

3.2.2 卡门涡街 随着雷诺数的进一步增高，圆柱后面的旋涡逐渐趋于非对称，在圆柱后面形成交错排列的卡门涡街。图5为Re=100时的流态演进示意图。可以清晰的看到上下两个旋涡生成、发展、脱落并向下游运动形成周期性交错排列的现象及尾流区流动的周期性摆动。
  
图3 Re=1、30时的方柱绕流流态



图5 Re=100时方柱绕流卡门涡街

  
图4 方柱绕流定常层流状态下回流区长度和雷诺数的关系 图6 方柱绕流Strouhal 数与Reynolds数的关系

由图可见从180.0 s开始至193.5 s为一个整周期。因此相应的St(fsU∞/d)数为0.148.本文计算了Re=60，80，100，120，150，180，200，250，300 的情况，并与Suzuki及Breuer的结果进行了比较，见图6.

3种计算结果虽然在数值上出入不大，另外在总体趋势上也是一致的。Re<150时St数随Re增大而增大；在Re=150时St数达到最大值；Re>150后St数随Re增大而减小，但当Re>250后St数基本上成为一个常数不再减小。

4 结论

建立了先进的并行计算的软硬件系统，实现了对浅水流动问题的并行计算方法研究。利用低Reynold数下二维方柱绕流算例对本文的并行计算程序进行验证，将流动形态、回流区长度及Strouhal数等特征参数与前人结果作了比较，结果表明本文提出的并行算法大大加快了收敛及求解速度，并且获得了良好的加速比和并行效率，可以进一步用于大型的天然河道流场的精确模拟。

在并行计算中，加速比和效率是并行程序进行评价的重要指标，其中，加速比定义为：单个节点上的运行时间和P个节点上运行时间的比，效率定义为加速比和计算节点个数之比。本文计算中，2台工作站运行时，加速比达到了1.79，效率为0.89；4台工作站运行时，加速比达到了3.35，效率达到了0.84.这证明本文并行程序有优秀的并行性，完全能够胜任大规模的数值模拟。

并行计算是解决大规模流场计算的必然途径，对我国防洪减灾、赤潮现象等问题的研究必将产生深远的影响。

参 考 文 献：

[1] Kashiyama K, Saitoh K , Behr M， Tezduyar T .Parallel finite element methods for large scale computation of storm surges and tidal flows[J]。Int .J。Numerical methods in fluids, 1997，(24):1371-1380 .

[2] MPI Forum. MPI: A Message Passing Interface[J]。International  Supercomputer Applications, 1994, 8(3/4): 165-414.

[3] Lawson C L, Hanson R J, Kincaid D R, Krogh F T. Basic linear Algebra Subprograms for Fortran Usage[J]。ACM Trans. Math. Softs,1979，(5): 308-323.

[4] Hammarling, McKenney A, Ostrouchov S, Sorensen D. LAPACK Users' Guide, 2nd edition[M], Philadelphia, PA: SIAM, 1995.

[5] Breuer M, Bernsdorf J, Zeiser T, Durst F. Accurate computations of the laminar flow past a square cylinder based on two different methods: lattice-Boltzmann and Finite-volume.

[6] Suzuki H, Inoue Y, Nishimura T， et al. Unsteady flow in a channel obstructed by a square rod (crisscross motion of vortex). Int.[J]。Heat Fluid Flow, 1993， 1(14): 2-9.

[7] George, P L.Automatic Mesh Generation: Application to Finite Element Methods[M]。U S A : WILEY, 1991.

[8] Miller G L， et al. Automatic mesh partitioning[A]。In: A. George, Gilbert J, and Liu J, eds. Sparse Matrix Computations: Graph Theory Issues and Algorithms. IMA Volumes in Mathematics and its Applications, Springer Verlag, 1993.

[9] 陈国良。并行计算[M]。北京：高等教育出版社,1999.

[10] Chun Bo Jiang，Qing Hai Zhang,Xiao Mi An. River flow simulations using parallel computing techniques[M]。XXIX IAHR Congress Proceeding ,Theme D.2001.