基于GPS的洪泽湖水流流场质点示踪及数值模拟

陈 钢， 翟 月， 王旭丹， 王 淳， 茅志兵

(1.河海大学 水文水资源学院， 江苏 南京 210098； 2.南京慧水软件科技有限公司，江苏 南京 210019； 3.天津市水文水资源管理中心， 天津 300061)

1 研究背景

湖泊作为静止水体具有保水、蓄水、清洁和供水的基本能力，是水文和生物地球化学水循环的重要组成部分，并影响着生态和经济等许多方面[1-2]。然而，众所周知，大型湖泊总是存在不同程度的水环境问题，比如富营养化严重、有机污染加剧、生物多样性遭到破坏等[3-5]。湖泊内水、沙及各种营养物质和污染物主要来源于流域，湖泊水流运动路径是各种入湖物质输移的主要途径与方式，并在湖内发生物理、化学和生物等一系列复杂反应，进而给湖泊水质的时空变化带来严重的影响[6]。

洪泽湖是我国第四大淡水湖，是淮河流域的最大调蓄湖泊及其下游地区重要的水源地，同时也是南水北调东线工程的过水通道，在防洪保安、城乡供水、交通航运以及维系生态平衡等方面有不可替代的作用[7-8]。但是，洪泽湖由于受粗放式发展模式的影响，生态环境退化问题愈发严重[9-10]。淮河上游及洪泽湖周边城市排放的工业污染物、生活污染物，洪泽湖区域内畜禽养殖与水产养殖过程中污染物直接入湖，不仅影响了洪泽湖综合功能的发挥，还使湖体健康受到严重威胁[11-15]。因此近年来，洪泽湖的健康和可持续发展问题引起了全社会的高度关注。水流流场与水质场的模拟是揭示湖水水流水质主要影响因素的有效工具，对湖泊的治理起重要作用[6]。尽管基于水动力学模型探究表层水流运动特性和空间分布特征的研究取得了一定进展，不少学者依据鄱阳湖等湖体构建了二维水动力数学模型来模拟风场条件下的水流时空结构[16-22]，但对于湖泊系统、从入湖口涌入湖区的水在空间尺度上的运动路径缺少实际验证。目前，只有李云良等[18]在鄱阳湖进行了野外粒子示踪实验，但由于投放粒子数量较少及监测时间较短，未能完全反映出湖泊整体的流场轨迹。

本文与先前研究的不同之处在于，基于GPS示踪实验和二维水动力模型，通过将野外实地监测与水流数值模拟相结合，深入探究了洪泽湖主要出入口及气象条件影响下的流场过程。本文在洪泽湖主要入湖口分别投放GPS示踪器，通过GPRS/CDMA网络实时传输，实现对表层水流轨迹、方向的实时跟踪，最终形成完整的洪泽湖初步流场图，同时构建二维水动力模型模拟验证基于GPS示踪器的流场的准确性。本文所进行的研究有利于正确认识湖泊上游来水的迁移路径、空间分布特征和湖区外来污染物的可能运动路径，为进一步深入研究洪泽湖水动力水质特征提供前提条件，为改善湖区生态健康提供可靠的科学依据，这些正是目前洪泽湖研究中极为重要又缺乏的方面。

2 材料与方法

2.1 研究区概况

洪泽湖位于江苏省西北部，淮安、宿迁两市境内，湖面辽阔，资源丰富，素有“日出斗金”的美誉，是淮河流域最大的湖泊型水库，地处苏北平原中部偏西，位于淮河中下游结合部[14]。湖泊最大长度约60.0 km，平均宽度约24.4 km，底部为浅碟型，总体呈西北高、东南低趋势。洪泽湖承泄淮河上中游15.8×104km2的来水[11]，入湖河流主要在湖西侧，主要有淮河、怀洪新河、池河、新汴河、老汴河、新(老)濉河等，其中淮河入湖水量占总入湖水量的70%以上。洪泽湖出湖主要包括淮河入江水道、入海水道、淮沭新河和苏北灌溉总渠，大部分的洪水由三河闸下泄入江水道后流入长江。洪泽湖是淮河流域最大的湖泊，由于湖面宽阔，在风场和上游入流的作用下，产生环流，对湖水流动有较大的影响，产生水面倾斜，影响大堤防洪安全和闸门的泄水流量。

本文的具体研究区域为洪泽湖湖区，其周边水系及主要出入口示意图如图1所示。洪泽湖主要的入流口有：淮河干流的蚌埠站、怀洪新河的峰山站、濉河的泗洪站、新汴河的团结闸等；主要的出流口有：三河闸、二河闸、高良涧进水闸、高良涧水电站等。这些出、入口都是基本的水文站，具有完整的水文资料，为水动力学模型的模拟计算提供了可靠的边界条件。

图1 洪泽湖周边水系及主要出入口示意图

2.2 GPS野外示踪实验

GPS以全天候、高精度、自动化、高效益等优点赢得广大测绘工作者的信赖，并成功应用于大地测量、工程测量、航空摄影测量、运载工具导航和管制、资源勘察、地球动力学等多个学科与领域。GPS的定位是通过围绕地球运转的人造卫星连续向地球表面发射连续波无线电信号，卫星导航接收机通过接受卫星发出的无线电信号，用导航算法解算得到定位仪器的位置[23]。本次示踪实验投放的GPS定位仪器由蓄电池、GPS用户接收机、天线等组成，可以每隔10 min向终端接收设备发送实时坐标。

野外实验前对淮河干流的蚌埠站、怀洪新河的峰山站、安河的金锁镇、濉河的泗洪站、新汴河的团结闸5大主要入口处及实际地形进行实地考察，最终选定4个主要入湖口处(图1)投放示踪器，同一入口的示踪器分为一组。为全面深入探究洪泽湖水体运动规律，针对不同入湖水文条件分别投放GPS示踪器。通过获取GPS示踪器标点的精确信息，结合可视化技术解析4组GPS示踪器在不同入湖水文条件下的实验结果，获得不同入口示踪器的真实运动轨迹，从而形成不同环境下的洪泽湖初步流场图。

2.3 二维水动力模型构建

水动力数学模型是将水动力学问题用数学方程进行描述[24]，并在一定定解条件下求解，从而模拟实际的水流状态[25-26]。目前数值模拟方法主要有：有限差分法(finite difference method, FDM)、有限体积法(finite volume method, FVM)和有限单元法(finite element method, FEM)。本文洪泽湖水动力模型采用二维浅水波方程来描述，洪泽湖是典型的二维区域，湖面宽阔，蓄水面积1 500 km2，因此必须考虑风应力和柯氏力的影响，采用完全的二维模型模拟，并用破开算子法分裂二维浅水波方程[27]，然后利用FDM方法求解：

(1)

(2)

式中：ρa为空气密度，kg/m3；CD为阻力系数；W为离水面10 m高处的风速矢量，m/s。

根据研究区域的资料情况，把流域概化为具体的数值计算模型。模型采用四边形网格(500 m×500 m)来划分单元，共计划分为6 174个计算单元(图2(a))，通过DEM(digital elevation model)数据生成湖底地形(图2(b))。本研究收集了2018年7月1日至2019年10月1日的洪泽湖出入湖流量水位及洪泽湖附近气象资料，设置模型的入湖流量边界为泗洪(老)、泗洪(濉)、双沟(合)、小柳巷、盱眙(图2(a))，设置出湖的水位或流量边界为二河闸、高良涧闸、三河闸(图2(a))。气象站资料采用距离三河闸2 km的蒋坝站(图2(a))，通过微型综合自动气象站连续监测降水、蒸发、风速、风向、气温等水文气象要素。根据湖泊水文资料设置的洪泽湖出入流边界条件信息如表1。

图2 洪泽湖水动力模型设置及湖底地形概化图

表1 模型边界设置情况

3 结果与分析

3.1 野外示踪实验结果与分析

本研究在2019年1月16日、2019年6月23日与2019年8月18日分3次投放示踪器，每次投放4组，每组投放5个示踪器。根据各组示踪器发送的位置信息解析4组GPS示踪器不同水文条件下的实验结果，观察时长共计约30 d。

通过对2019年1-9月洪泽湖流量边界的分析(图3)得出，洪泽湖主要的入流为盱眙和小柳巷，因此可视为控制条件。2019年1月16日第1次投放后的30 d，小柳巷与盱眙的平均入流量分别为312 m3/s和815 m3/s；2019年6月23日第2次投放后的30 d，小柳巷与盱眙的平均入流量分别为560 m3/s和1 150 m3/s；2019年8月18日第3次投放后的30 d，小柳巷与盱眙的平均入流量分别为131 m3/s和580 m3/s。因此在本研究中，通过这3次投放周期的平均流量比较，假定第1次投放为中流量入湖时期，第2次投放为大流量入湖时期，第3次投放为小流量入湖时期。

图3 2019年洪泽湖入湖流量及示踪器投放时间

3.1.1 第1次投放后水流轨迹解析 2019年1月16日进行示踪实验的第一次投放，图4展示了中流量入湖条件下示踪实验仪器自投放后5、10、20及30 d的轨迹，其中各示踪器运动轨迹使用蓝色线条绘出，当日示踪器具体位置使用红点标明。第1组的示踪器在为期30 d的示踪过程中没有发生位置的变动，即湖泊西北凹岸处入流的水流处于静止状态，没有向湖泊中心运动。第2组的示踪器轨迹显示，自湖泊西北部入流的水流首先向南运动，10 d后到达双沟镇附近并转为向东运动，20 d后到达临淮镇，此后水流向东北运动，在投放30 d后到达龙集镇南部附近，整体运动较为缓慢，预计此后水流向东运动。第3组的示踪器运动轨迹显示，自湖泊西南部双沟镇附近入流的水向东北方向运动，5 d后到达临淮镇附近，并在投放后20 d到达湖泊中心，30 d后从湖泊东部出流口的二河闸流出，运动轨迹受湖泊北岸岸线形状的影响，呈现为向北凸出的曲线。第4组示踪仪器的轨迹则表明湖泊南部入流口的水流向东北流动，此后在湖泊南部凸岸处变为东向运动，其迹线随时间推移而愈发分散，在投放后20 d已有水流自湖泊东南部出口的三河闸流出，此时水流在湖泊东南岸的东双沟镇附近分为两部分，小股水流向西北方向运动，在30 d后从东部二河闸流出，水流主体向南流动，在30 d后全部从湖泊东南部的三河闸出湖。

图4 示踪实验仪器第1次投放后的运动轨迹

3.1.2 第2次投放后水流轨迹解析 2019年6月23日进行了第2次投放，图5展示了大流量入湖条件下示踪实验仪器自投放后5、10、20及30 d的轨迹，其中各示踪器运动轨迹使用黑色线条绘出，当日示踪器具体位置使用红点标明。由编号为第1组的示踪器轨迹可得，自湖泊西北凹岸处入流的水流缓慢地向南运动，相比中流量入湖条件的迹线水流状态由静止变为运动，可以推断大流量条件下，湖泊西北部的局部流速比中流量条件的流速更大。第2组示踪器30 d的轨迹显示自湖泊西北部入流的水流向南运动，相比中流量条件下的轨迹基本相似却略显动力不足。7月22日，即运动30 d后仅到达临淮镇，预计此后水流继续向东运动。结合第1组的轨迹可以推断，湖泊西北部的局部流场较为复杂，流速的分布与中流量条件下的情况相差较大。第3组的示踪器轨迹则可以看出，自西南部双沟镇附近入流的水向东运动，此后迹线向北凸的程度明显降低，而更趋近于横跨湖泊中心的直线，10 d后到达湖泊中心区域，在投放后20 d从湖泊东部出流口的二河闸流出。第4组示踪仪器的轨迹则表明湖泊南部入流口的水流向东北方向流动，在投放后5 d即到达中流量条件下该组仪器投放20 d后到达的位置，并已发生水流分向，分为西北向流动的水流主体与南向流动的小股水流，西北流向的水流在投放后10 d即到达湖泊东部出流口的二河闸，向南流的小股水流也在10 d后从湖泊东南部的三河闸出湖，其运动轨迹相比第1种流量条件更为分散且弯曲程度更大。

根据与图5呈现的中流量入湖条件下水流运动轨迹的比较可以发现，大流量条件下湖泊西北部的局部水流运动更为复杂，而南部水流运动明显加快，整体水流从入流到达湖泊出流口的时间明显缩短；另一方面，整体的水流运动轨迹与中流量条件下的水流轨迹相比有一定程度的改变。自湖泊西北部入流的水流流动轨迹基本与中流量条件一致；自西南部入流的水在横跨湖泊中心区域后，迹线向北凸的程度明显降低，而更趋近于横跨湖泊中心的直线；湖泊南部入流口水流的迹线更为分散，主体水流由向南运动变为向北运动。

图5 示踪实验仪器第2次投放后的运动轨迹图

3.1.3 第3次投放后水流轨迹解析 2019年8月18日进行第3次投放，图6显示了该小流量入湖情况下，示踪实验仪器自投放后5、10、20及30 d的轨迹，其中各示踪器运动轨迹使用橙色线条绘出，当日示踪器具体位置使用红点标明。第1组示踪器在静止了10 d后开始向南运动，30 d后的运动距离相比大流量条件稍短。第2组示踪器的运动轨迹与前述两种流量条件下的轨迹相似，30 d后到达的位置则在湖泊北岸的半城镇附近，处于相同时间中流量与大流量条件下到达的位置之间，预计此后将继续向东运动。第3组投放的示踪器轨迹则显示，自湖泊西南岸双沟镇附近入流的水向东流动，10 d后到达临淮镇附近并变为向东北流动，在20 d到达龙集镇附近的湖泊中心后继续向东运动，此后运动轨迹为向北凸出的曲线，30 d后到达湖泊东岸的出流口二河闸，运动轨迹与第1种流量条件下的运动轨迹一致。第4组示踪仪器的轨迹表明湖泊南部入流口的水流向东北运动，8月27日到达老子山镇的东北侧附近并改变流向为向东流，此后与前两种流量条件不同，没有发生水流分向现象，到达东双沟镇附近后向南运动，在30 d后自湖泊东南部出口的三河闸流出，其余示踪器受动力不足的影响部分在入流口附近搁浅，或仅到达老子山镇北侧附近的位置。

图6 示踪实验仪器第3次投放后5d、10d、20d及30d运动轨迹图

根据第1～3组的示踪器轨迹，可以看出湖泊中心流速较大，而洪泽湖南部绕盱眙县北部一段的流速较前两种流量条件明显降低。小流量条件下水流运动会受到湖泊表面风场的一定影响，洪泽湖南部的水流在9月6日后持续受到北风的影响，向南的水流运动轨迹不再是曲线而是接近直线。

3.1.4 3次示踪实验综合分析 与中、大两种不同流量入湖条件下的水流轨迹进行比较，小流量入湖条件下湖泊西北部的水流运动情况十分复杂，流速时空差异较大，而南部的水流运动则随着流量变小明显变得迟缓；另一方面，整体的水流运动轨迹与中流量入湖条件下的水流轨迹更为相似。

总结上述推断可以发现，流经湖泊中心区域的水流运动轨迹均为自西向东的曲线，均从二河闸流出，即湖泊主要的水流运动是自西向东，主要出口是二河闸。而自湖泊西北角以及南岸入湖的水受到流量条件的影响，流量较小时水流难以进入湖泊的中心流场，因而流动缓慢或处于静止状态；流量较大时流速空间分布不均，流场较复杂。例如南岸入流的水流有两个出流点，流量大时受湖泊主流向东运动的影响，发生水流分向，流量越大，参与湖泊水体整体向东运动的水越多，在示踪实验结果中体现为监测到从二河闸流出的示踪器数量越多，从三河闸流出的示踪器数量越少；流量越小时，水流受湖泊南岸附近局部流场的影响越大，水流流动缓慢或停滞，仅小股从湖泊东南岸的三河闸流出。这也与三河闸不常开启，湖泊水流较少从该处流出的现象相符合。

此外，小流量条件下水流运动会受到湖泊表面风场的一定影响，例如洪泽湖南部的水流在9月6日后持续受到北风的影响，向南的水流运动轨迹不再是曲线而是接近直线。而流量较大时水流流动的主要影响因素是湖泊自身的流场，湖泊表面风场的作用较小。

3.2 水动力模拟结果分析

采用二维水动力模型模拟洪泽湖二维水流动态。该模型使用2018年7月1日至2019年10月1日出流口二河闸、三河闸的水位与流量资料进行率定与验证。站点拟合的确定性系数可达0.91，相对误差均在8%以内，表明模型能够很好地模拟洪泽湖内水位的时空变化(图7)。图7(b)显示了2019年7月2日洪泽湖表面的二维流速场，其中代表流速的箭头长度越大表示流速的数值越大。由图7(b)可以明显看出，湖泊的西部(靠近湖泊西部入流口的临淮镇)、湖泊南岸、湖泊东南部入江水道及东部二河闸两个出口处的流速较大，其中湖泊东南部的二河闸的流速最大，这与二河闸作为洪泽湖的主要出流点的实际情况吻合。湖泊中心区域流速较大，整体呈现自西向东的流动特征，也与前述示踪实验的结果基本一致。

图7 洪泽湖流速场模拟结果与示踪迹线叠加

为综合本研究质点示踪实验与水动力模拟这两种研究手段的优点，并更直观地展示洪泽湖水流流场的特性，将上述中、大两种流量条件下示踪器实验投放10 d后的结果与水动力模型模拟的当日流场分别叠加，得到图7(a)、7(b)所示的结果。

由图7可见，中流量入湖条件下，湖泊南岸流速较大，因此南岸入流的水流开始时主要是向北运动，随后受到湖泊中心流场的影响而转为向东流动，唯有一股水流由于靠南岸较近，未参与中心流场的水流运动，短暂地向东流动后迅速在老子山镇附近转为向南流动。而湖泊西北角的狭长水域，水流流速差异较大但总体较小，因此同属于西北岸入流，虽然入流点位置不同的水流其运动轨迹相似，但运动距离不同。湖泊西岸临淮镇附近水流流速大，是湖泊主体水流向东运动的起点，从此处入流的水运动速度快，并且流经湖泊中心区域直接到达湖泊东岸的二河闸出口。

大流量入湖条件下，中心流场的范围扩大，流速整体增大，这体现在示踪实验中即湖泊各处入流的水流运动均明显加快，除小股水流受南岸局部复杂流场影响最终从三河闸流出外，参与湖泊主体向东流动的水流在10 d内全部到达二河闸出口，湖泊西北部的水流运动也有整体加速的趋势。

小流量入湖条件下，中心流场范围缩小，整体流速降低明显，西北部流速空间分布不均，而湖泊南岸附近流速明显变小。示踪实验中自西北部不同地点入流的两股水流运动轨迹完全不同，一股水流静止而另一股水流运动速度较快，正是水动力模拟结果西北部流场复杂情况的体现。而湖泊南岸由于流速明显变小，示踪实验中水流运动明显变缓，相比前述的两种流量条件，其运动轨迹更靠近南岸，因此受南岸局部流场的影响较大。

综合上述3种情况可以看出，不同流量大小的入湖条件下，湖泊主体水流向东流动的速度随之变化，流场范围也与流量大小呈现正相关关系。除此之外，湖泊西北部的狭长区域与湖泊南岸附近的流场较为复杂，局部流速时空差异较大，与示踪实验所示的水流轨迹较为混乱的特征基本符合。水动力模型模拟结果的流场基本符合水流质点GPS示踪实验的结果，两者拟合较好，模拟较为可信。

4 结论与讨论

为探究洪泽湖水流运动及其具体流场分布，本文分别进行了以GPS为基础的水流轨迹示踪实验与二维水流数值模拟，主要结论如下：

(1)示踪实验中，流经湖泊中心区域的水流运动轨迹均为自西向东的曲线，均从二河闸流出，可以推断，湖泊主要的水流运动是自西向东，主要出口是二河闸。湖泊西北角以及南岸入流的水受到流量条件的影响，流量较小时水流难以进入湖泊的中心流场；流量较大时流速空间分布不均，流场较复杂。南岸入流的水流有两个出流点，流量大时受湖泊主流向东运动的影响，发生水流分向，流量越大，参与湖泊水体整体向东运动的水越多；流量越小时，水流受湖泊南岸附近局部流场的影响越大，水流流动缓慢或停滞，仅小股从湖泊东南岸的三河闸流出。此外，小流量条件下湖泊整体水流运动会受到湖泊表面风场的影响，而流量较大时水流流动的主要影响因素是湖泊自身的流场，湖泊表面风场的作用较小。

(2)将不同流量条件下的示踪实验结果与水动力模型模拟结果相叠加，不同流量条件下，湖泊主体水流的向东流动速度随之变化，中心流场范围也与流量大小呈现正相关关系。除此之外，湖泊西北部的狭长区域与湖泊南岸附近的流场较为复杂，局部流速时空差异较大，与示踪实验所示的水流轨迹较为混乱的特征基本符合。水动力模型模拟结果的流场基本符合水流质点GPS示踪实验的结果，两者拟合较好，数值模拟较为可信。