基于Google Earth的分蓄洪区水沙模拟与演示系统研究

，

(长江科学院 河流研究所，武汉 430010)

自然溃堤或人工分洪情况下，洪水在洪泛区内的演进将给洪泛区内人民生命财产安全带来巨大威胁，因此模拟洪水在洪泛区内的演进过程，将给防洪决策和防洪抢险提供技术支撑。目前，许多学者对分蓄洪区洪水演进开展了数值模拟研究[1-3]，然而分蓄洪区地形获取的不便和模拟成果显示不直观等问题，仍然是制约数学模型发展的重要因素。近10多年来，地理信息技术得到快速发展，形成了许多的地理信息平台。因此，如何将洪水演进模拟结果显示在地理信息平台之上，是一件非常有意义的事情。

Google Earth软件(以下简称GE)便是其中应用最为广泛的通用的地理信息平台之一，它把卫星照片或航空照相和GIS信息布置在一个三维的地球模型上，而且免费对大众开放。本文基于Google Earth软件平台，建立了分蓄洪区水沙演进数学模型，能快速完成模型数据前处理和模拟计算，并实现了计算结果在三维虚拟场景中的动态演示。该模型系统适用于紧急条件下的快速计算和直观展示，计算成果可为分蓄洪区防洪决策与防洪抢险提供技术参考。

1 基于非结构网格的数值求解方法

分蓄洪区一般不规则，为此采用非结构网格[4]进行分蓄洪区的网格划分是非常合适的。非结构网格可以是三角形网格、四边形网格或三角形与四边形混合网格。

模型基本方程采用守恒形式的水深平均平面二维水流基本方程[5]。采用基于同位网格的有限体积算法进行方程的离散[6]，流速矢量基于直角坐标系，节点布置采用单元中心法；采用SIMPLEC算法进行离散方程的求解，该算法的优点在于即便采用较粗的网格也能保持水量和动量的守恒。

2 基于Google Earth的模拟演示系统

本文基于Google Earth软件进行二次开发，将数值模拟成果在由航拍地貌和三维建筑物构成的虚拟地球上进行显示，有助于更好地分析计算成果的合理性，同时还可利用Google Earth软件缩放、平移、鸟瞰等功能，更直观地查看计算成果。

Google Earth软件的二次开发有2种方式，一种是基于ComAPI；一种是基于KML(全称为Keyhole Markup Language)。基于 ComAPI的开发方式主要用来控制GE的视角,实现动画效果，而基于KML的开发方式主要用来生成地理要素(如点、线、面、多边形、多面体以及模型等)，实现数据的动态更新等。在实际应用中，需要结合这2种开发方式，才能得到比较理想的效果。

笔者采用VB和Fortran语言混合编程方法，对Google Earth软件进行二次开发，研制了“基于虚拟地球的水沙水质数值模拟与演示系统V2.0”，系统由前处理、模型库、后处理及数据库等部分构成，模型可实现干河床的水沙演进计算。该系统不仅实现了界面操作的可视化，而且实现了计算结果在虚拟地球上的动态演示，系统主界面见图1。

图1 可视化界面图

2.1 基于ComAPI的开发

GE 公布的API目前是1.0a版，数量比较少，所能实现的功能也不多，在实际开发中主要用来加载KML数据、控制视角和实现动画播放。加载KML数据命令有：

(1) IApplicationGE. OpenKmlFile(kmlFilePath)，加载指定的KML文件到GE中；

(2) IApplicationGE.LoadKmlData(kmlData)，加载KML文本流到GE中。

2.2 基于KML文件的开发

2.2.1 坐标转换

目前我国的河道地形图大地基准面一般采用54坐标系或西安80坐标系，而Google Earth软件是基于WGS84坐标系统，因此在生成KML文件之前需对计算结果文件进行坐标转换。上述3类坐标系的转换在相关教科书中有详细的说明，这里不再赘述。

2.2.2 流速矢量场

计算成果显示包括流速矢量场和水位、水深、流速等物理量的标量场，同时还需文字标注、图例。为进行流场的动态演示，需进行每帧流场的时间控制。

图2 流速箭头示意图

计算域内每个点的流速矢量可用箭头表示，箭头长短表示流速大小，箭头方向代表流速方向。每个箭头可由6个有序点形成的折线表示，如图2所示。

在KML文件中，每个点的流速矢量可用LineString标签来定义，代码如下：

112.138 084 0, 30.284 331 1, 0

112.140 711 6, 30.283 183 1, 0

......

其中字段表示定义标签；定义该箭头名称；定义流速箭头颜色，需根据流速的大小预先定义；绘制折线；下面定义6个点的坐标，注意这里采用的是经转换后的经纬度坐标。

2.2.3 等值线云图

等值线的绘制与流速箭头的绘制类似，只需绘制每条折线就行，云图的绘制需将封闭折线构成的区域进行颜色填充，可用Polygon标签来定义，代码格式如下：

112.141 431 5, 30.283 587 2, 0

112.140 788 8, 30.281 146 2, 0

112.144 288 4, 30.282 297 2, 0

112.141 431 5, 30.283 587 2, 0

112.141 431 5, 30.283 587 2, 0

2.2.4 文字标注

流速值、水深值或其他文本的标注可采用Google Earth中的Point地标标签来定义，代码如下：

2.2.5 流场图例

完整的流场显示需定义图例，在Google Earth中图例由ScreenOverlay标签控制，在动态演示过程中保持不变，代码如下：

2.2.6 时间控制

为实现流场的动态演示，需定义每帧流场的起止时间，采用TimeSpan标签定义，代码如下：

3 在荆江分蓄洪区的应用

荆江分蓄洪区位于长江南岸，是长江中游防洪系统的主要组成部分，分洪区面积921 km2，有效蓄洪容量54亿m3。荆江分洪工程包括进洪闸、节制闸和分洪区围堤。进洪闸位于太平口东岸，设计进洪流量8 000 m3/s。荆江分蓄洪区的功能在于分泄荆江河段的超额洪水，使洪水能够安全通过，以保护荆江两岸广大农田和人民生命财产安全。

1954年7月22日至8月1日荆江分洪区进行了3次分洪运用，分洪总量为122.6亿m3。采用建立的模型系统，对1954年荆江分蓄洪区分洪过程进行了模拟计算。

3.1 计算条件

计算范围为整个荆江分蓄洪区，以北闸为进流边界，南闸为出流边界。网格采用非结构网格，见图3。计算网格地形根据航拍的DEM地形数据插值得到。

图3 分洪区计算网格布置

将荆江分蓄洪区分洪过程概化为溃坝过程，溃坝计算关键在于如何进行干河床与动边界的处理，笔者提出了“干湿双最小水深法[6]”，能较好反映溃坝水流在干河床上的演进过程。

由于荆江分蓄洪区内缺乏实测的洪水资料，因此，糙率参考相关资料确定[6]：树林0.070，旱地0.065，水田0.050，水面0.025。如果某网格内含有多种地形，则按照各种地形糙率的加权平均值确定该网格的糙率。

模拟计算从1954年7月22日8:00开始，考虑为全溃，直至蓄满为止。进口含沙量假定为零，床沙中值粒径取0.2 mm。

3.2 计算成果分析

3.2.1 洪水演进计算成果

计算结果表明，分洪闸下泄流量先缓慢增大，最大可达8 000 m3/s，后快速减小，直至蓄满所需时间为11 d左右。

图4给出不同时刻洪水在分蓄洪区的演进情况(图中箭头表示流速矢量大小和方向，蓝色云图表示水深，蓝色越浑则水深越大)。从图4中可见，系统较直观地反映了洪水在航拍地貌上的演进过程，在溃口和分洪区缩窄段，流速明显增大，最大流速为1.2 m/s，图中彩色云图能直观反映不同时刻洪水淹没范围及淹没水深。

图4 不同时刻洪水演进流场图

3.2.2 泥沙输移计算成果

图5(a)给出了分洪67 h时分蓄洪区含沙量分布情况(图中颜色越红表示含沙量越大)。从图4中可见，由于进口未给定含沙量，泥沙主要来自于床面的冲刷，因此流速越大部位(如图5(a)中深蓝色区域)，含沙量也会越大，最大含沙量可达2 kg/m3左右，位于分洪口门下游附近；而在分洪区下游，由于流速不断减小，含沙量也随之减小。

图5 分蓄洪区含沙量和冲淤厚度分布(T=67 h)

图5(b)给出了分洪67 h后分蓄洪区冲淤厚度分布情况(图中红色表示淤积，蓝色表示冲刷，颜色越深则冲淤幅度越大)。从图中可见，一般大流速部位，床面易发生冲刷，最大冲深为1.7 m，其余大部分区域发生淤积，最大淤厚约1.5 m。

4 结 论

(1) 建立了基于地理信息技术的分蓄洪区平面二维水沙演进数学模型。模型能直接在Google Earth软件环境中确定计算域边界，采用无结构的三角形、四边形网格进行计算域的网格划分，根据航拍的DEM地形插值得到计算网格地形高程，提高了前处理的效率。

(2) 对Google Earth软件进行二次开发，研制了计算结果后处理模块，实现了分蓄洪区水沙演进过程在三维航拍地貌中的动态演示，并能与Google Earth同步实现缩放、鸟瞰、漫游等多种查看功能。

(3) 该模型系统在长江荆江分洪区得到成功运用，复演了1954年分洪后的水沙演进过程，并在Google Earth地貌中进行动态演示，可直观查看水深、流速、洪水淹没范围、冲深等要素，可为洪水风险评估和应急预案的制定提供参考。

参考文献：

[1] MACCHIONE F, MORELLI M. Practical Aspects in Comparing Shock-Capturing Schemes for Dam-break Problems [J]. Journal of Hydraulic Engineering, 2002, 129(3)：187-195.

[2] 余明辉，张小峰. 平面二维溃堤水流泥沙数值模拟[J]. 水科学进展，2001， 12(3)：286-290.(YU Ming-hui, ZHANG Xiao-feng. Horizontal 2-D Uneven Sediment Mathematical Model in Dike Burst[J]. Advances in Water Science, 2001, 12(3): 286-290.(in Chinese))

[3] 王志刚，汪德爟，赖锡军, 等. 下游为干河床的溃坝水流数值模拟[J]. 水利水运工程学报，2003，(2):18-23.(WANG Zhi-gang, WANG De-huan, LAI Xi-jun,etal. Numerical Simulation of Flow on Dry Bed for Dam Break[J]. Hydro-science and Engineering, 2003, (2):18-23.(in Chinese))

[4] 张细兵.河道有限元网格自动剖分方法研究[J]. 长江科学院院报，2002，19(3):19-21.(ZHANG Xi-bing. Auto-generating Method of Triangulation Grids for River Reach[J]. Journal of Yangtze River Scientific Research Institute, 2002,19(3):19-21.(in Chinese))

[5] 张细兵，范北林. 溃坝洪水演进平面二维数学模型初步研究[J]. 长江科学院院报, 2006，23(6):19-21.(ZHANG Xi-bing, FAN Bei-lin. Numerical Simulation for Two-Dimensional Dam-break Flood Flow[J]. Journal of Yangtze River Scientific Research Institute, 2006, 23(6): 19-21.(in Chinese))

[6] 张细兵，欧治华，崔占峰，等. 基于非结构网格的分蓄洪区水沙演进数学模型研究[J].长江科学院院报，2011，28(4):75-79.(ZHANG Xi-bing, OU Zhi-hua, CUI Zhan-feng,etal. Unstructured Grid Model of Flow and Sediment Evolution in Flood Diversion and Storage Area[J]. Journal of Yangtze River Scientific Research Institute, 2011, 28(4):75-79.(in Chinese))