开放科学(资源服务)标识码(OSID): 
0 引言
近年来,随着数字孪生技术的快速发展及智慧水利建设的深入实施,水利部启动了数字孪生流域先行先试工作,在大江大河重点河段、重要水利工程、主要支流等区域开展数字孪生建设[1-2]。水质安全的预警与管理是数字孪生流域建设的重要环节,而对水污染源的扩展过程进行准确数值计算是开展模拟推演的基础[3-4]。长期以来,国内外研究机构对各类水域的水动力过程及水质变化机理开展了富有成效的研究,形成了一大批考虑因素全面、功能完善的环境水动力与水质耦合模型和计算软件,其中最具代表性的有:美国国家环境保护局支持开发的环境流体动力学(Environmental Fluid Dynamics Code,EFDC)模型、荷兰Delft水力学实验室开发的Delft3D模型软件、丹麦水资源与环境研究所(Danish Hydraulic Institute,DHI)开发的Mike系列软件等。这些模型软件可以模拟河流、湖泊、水库、湿地、近岸海域等区域的泥沙传输、污染物扩散、水体富营养化等过程,目前已在高校、政府和工程及环境咨询机构中广泛使用[5-6]。
利用各类水质模型软件计算出水质过程后,为直观呈现和分析水质演化趋势,需要借助三维可视化工具对计算结果进行可视化。部分水质模型计算的商业软件自带简易的三维可视化模块,可以在计算出水质结果后直接进行可视化,如EFDC Explorer、Delft3D等,但其可视化功能单一,缺乏统一的空间坐标基准和数据转换集成能力,因而无法接入地理空间、水文测站、天气环境等数据,仅能通过网格分级设色的方式呈现水质模型本身,难以对水质演化过程形成系统、全面的认识。
为了实现水质计算结果与其他各类信息的集成,以地理信息系统(Geographic Information System,GIS)技术为基础的开发框架被广泛使用,研究者利用GIS引擎的空间数据融合能力进行扩展开发,对水质模型结果进行转换与集成。GIS引擎借助空间数据分级、索引构建、分页调度等机制实现静态数据的快速集成和调度,但针对动态数据的图形渲染和更新能力较弱、显示效果不佳[7],因而基于GIS引擎的水动力-水质可视化研究大多采用二维平面展示方式,针对水质分层中的某一层级、极值或均值指标等进行分级设色渲染,呈现各网格的水质分布情况,并采用时间轴播放的方式模拟水质演化过程[8-9]。杨坤等[10]基于MapTalks地图引擎建立了乌东德水库水质模拟推演系统,模拟化学需氧量、总磷、总氮等水质因子的沿程变化,采用二维地图动态展示各网格单元的水质状况;邹友琴等[11]基于桌面端GIS组件设计了面向鄱阳湖水质数据模拟和分析的软件工具,开展有机污染物的扩散模拟和趋势预测;吴豪杰[12]基于ArcGIS Engine搭建了清潩河水质三维仿真系统,利用序列快照播放方式对水体表层和重要断面的污染物浓度进行动态渲染模拟。
随着可视化仿真技术的快速发展和数字孪生引擎的升级迭代,研究者尝试基于三维空间插值、有限元分析、点云阵列拟合、体素模型绘制等新技术开展三维实体数据的渲染研究,为三维水质推演仿真提供了新的思路。赵纪堂等[13]基于改进的反距离权重(Inverse Distance Weighted,IDW)插值法构建了淮南市顾北矿新生界下部含水层的真三维地下水化学场,为水质真三维空间分布规律研究提供辅助决策支持;徐劲力等[14]基于Unity3D引擎开发了有限元分析结果的虚拟评审应用,实现基于混合现实场景的有限元分析结果的可视化,为三维水质模拟、模型受力分析、变形仿真等类似应用场景的计算模拟提供了可视化思路;刘振华[15]提出了基于水质变化数据的时空体素模型构建方法,并运用光线投射算法实现高效体绘制从而模拟水质变化。这类研究主要侧重于可视化技术本身,尽管目前缺少在三维水质推演仿真领域的相关应用,但其多样化的仿真方法和手段亦可为水质仿真提供另类解决方案。
现有的可视化方案能在一定程度上直观展示水质模型计算结果,但存在下列问题:
(1)目前常用的水质分析专业软件仅提供简易的三维可视化模块,只能对水质网格模型本身进行可视化,无法接入高精度地形、实景三维模型等,难以实现水文气象、防洪形势、工程调度等数据与数字孪生场景的集成融合。
(2)现有的水质推演研究主要以GIS引擎为基础,受制于引擎的渲染性能和仿真效果,可视化方案大多为二维平面展示方式,局限于展示水质分层中的某一层级、极值或均值指标,难以表达整体的水质状况。
(3)部分可视化研究通过网格关联映射、点云阵列拟合、时空体素模型等方法实现水质模型的三维可视化,通常将水体区域采样为离散的点云阵列,或划分为相互独立的垂向立方体,相邻区域的连接不平滑,导致三维展示结果存在立方体的颗粒感,且颜色分级渲染不连续,可视化效果不佳。
鉴于此,本文借助数字孪生技术对物理世界表达的优势,研究三维水质模型计算结果与数字孪生场景的集成模式,提出面向数字孪生场景的三维水质推演仿真方法,实现三维水质推演过程的直观、准确、流畅、动态呈现。
1 三维水质推演仿真方法
1.1 总体技术路线
为实现三维水质模型计算结果在数字孪生场景中的推演仿真,本文方法的总体技术路线如图1所示。该技术方案围绕水质推演的数据底板和业务应用需求搭建数字孪生场景,并在数字孪生场景下开展水质计算边界条件设定、计算模块调用、计算结果优化处理、结果集成和推演仿真、孪生场景中的调度方案/工程应用调整、计算边界条件更新迭代的全过程应用。
数字孪生场景构建是开展水质推演的基础,在完成数字孪生仿真引擎的选型和基础地理空间数据的集成后,需进一步结合业务需求集成水文气象、防洪形势、工程调度等数据构建数据底板,并实现模型边界条件设定、仿真推演控制、场景信息联动等动态交互操作功能。
水质计算结果的推演仿真方法是本文探讨的重点,共包含三大步骤:三维水质模型计算结果提取解析、推演网格与属性信息组织优化、数字孪生网格体构建与推演仿真。其中“三维水质模型计算结果提取解析”步骤将计算结果的二进制文件解析转换为水质推演结构化时序数据,为数字孪生仿真场景提供三维水质数据输入;“推演网格与属性信息组织优化”步骤从水质网格平面和垂向坐标优化、水质网格时序属性信息优化水质网格边缘矢量线提取等方面进一步对数据进行优化组织,降低数据冗余、提升展示效率和效果;“数字孪生网格体构建与推演仿真”步骤在数字孪生仿真引擎中构建三维水质孪生网格体,将属性信息与之关联映射,依据污染物浓度按照分级色带平滑渲染,提供多种可视化模式,并按照时间序列进行推演仿真。
1.2 三维水质推演数字孪生场景构建
为构建水质推演的数字孪生场景,需首先开展数字孪生仿真引擎的选型。鉴于水质推演业务既需要场景效果逼真,能直观动态呈现水质污染过程,又要求数据和结果准确真实,能通过统一的地理空间坐标实现物理世界向虚拟世界的映射,因而本研究选用Unreal Engine 5游戏引擎+Cesium For Unreal插件作为数字孪生仿真引擎的基础。其中Unreal Engine 5提供三维底层数据结构、环境光照、物理引擎、渲染特效等支撑能力,Cesium For Unreal插件提供全球尺度的地理空间场景数据接入、坐标转换与多尺度数据调度能力。
完成数字孪生仿真引擎的选型并进行功能扩展开发后,需要进一步构建水质推演的数据底板,集成地理空间、水文气象、防洪形势、工程调度等数据,并添加用于态势与机理分析的模型库和规则库。在数据底板之上,水质推演的数字孪生场景还需提供场景信息交互功能。一方面可在场景中交互查询和设置水质模型计算所需的污染物类型、空间位置、计算时段、计算范围等边界条件,并将参数传递给水质计算模块;另一方面可集成水质推演仿真结果并以时间轴驱动场景信息的更新,进而在推演过程中交互式调整工程设施运用等条件,从而更新水质计算的边界条件并重新迭代计算。
1.3 三维水质模型计算结果提取解析
基于选取的边界条件计算得到水质推演结果的二进制文件后,需从中解析水质网格空间信息与水质时序属性信息。表1以EFDC模型为例列出了需要解析的文件。
表1 三维水质模型计算结果文件Table 1 Computed result files of the 3D water quality model |
| 文件 | 记录内容 | 详情 |
|---|---|---|
| Corners.inp | 网格坐标 | 水质正交网格顶点平面坐标 |
| EE_WC.OUT | 河底高程 | 各网格水下地形的河底高程 |
| EE_WS.OUT | 分层水深 | 各网格时序分层水深 |
| EE_VEL.OUT | 流速流向 | 各网格时序分层流速、流向 |
| EE_WQ.OUT | 水质指标 | 各网格时序分层污染物浓度 |
水质推演计算结果中包含网格的空间信息和时序属性信息。对于平面坐标、河底高程等空间信息,将其统一转换为WGS1984坐标,与数字孪生场景统一坐标基准;对于分层水深、分层流速、分层流向、分层污染物浓度等时序属性信息,将其按照水质推演时刻、水质网格所在层级、网格编号组织为三维数组dataArr,水质推演结果帧数为T,水质网格垂向层数为L,则时刻t分层l网格编号n的数据为dataArr[t][l][n]。
1.4 推演网格与属性信息组织优化
1.4.1 水质网格坐标优化
三维水质通常采用正交四边形网格建模,各水质网格由4个顶点组成,记录了网格的x和y方向投影坐标,对所有网格的顶点进行共用计算,共用顶点的 与 的容差应满足下列条件,即
式中: 分别表示顶点 在x和y方向上的空间投影坐标; 分别表示顶点 在x和y方向上的空间投影坐标。
若满足上述条件则认为两顶点相同,记录该顶点涉及的网格编号。对水质推演网格中相互邻接的多边形分别进行处理,消除因按空间单元划分产生的水质网格在接边处的高程值和信息差异,形成优化后的平滑连续网格信息如图2所示,并将其记录为水质网格顶点和拓扑连接数据。
1.4.2 时序属性信息组织优化
考虑三维水质的相邻网格划分大小基本一致,可认为每个网格的属性信息对顶点的权重相等,则对于某共用顶点V,该点位置的总磷浓度计算值为V(TP)。依据计算得出的顶点共用关系,其涉及的网格数量为N,网格编号分别为 ,其中 ,以总磷水质属性TP为例,则优化后的顶点水质信息为
式中Gi(TP)表示网格Gi处的总磷浓度。
按照此方式对水位、流速、流向、分层水深、其他污染物分层指标属性进行优化计算。
1.4.3 水质网格边缘矢量线提取
在三维水质推演过程中,由于内部水质网格被遮挡,若独立绘制所有水质网格边缘,会极大降低渲染和显示效率。因此可根据网格间的拓扑关系,提取水质网格的边缘矢量线并构建封闭边缘面,在保证显示效果的同时提高渲染效率。
经顶点共用计算后,获得了每个共用顶点涉及的网格数量。根据三维水质模型构网的规则,若顶点涉及1~2个网格则为边缘点,若涉及4个网格则为内部点,若涉及3个网格可能为边缘点或内部点。分别对涉及1、2、3个网格的顶点进行标记,并基于这些特征点遍历计算获取边缘,其过程如图3所示。
定义只涉及单个网格的顶点为角点,以其为特征点进行遍历计算,包括下列步骤:
(1)定义矢量线的顶点序列为 ,选取某顶点V,其涉及的网格数N(V)=1,以该点为起点进行遍历,并将点V加入到 顶点序列中。
(2)顶点V所在的网格为G1,获取该网格4个顶点的连接顺序,则顶点在该网格中的左右邻接顶点分别为VL和VR,获取它们涉及的网格数N(VL)、N(VR),若N(VL)≤N(VR)且N(VL)<4,则将点VL加入 顶点序列中。
(3)若N(VL)=2,该点为边缘点,获取VL涉及的网格G1、G2,在G2网格中,顶点VL在该网格中的左右邻接顶点分别为VLL和VLR,按步骤(2)的逻辑继续遍历;若N(VL)=3,为边缘点或内部点,获取VL涉及的网格G1、G2、G3,以G2和G3网格为计算单元,按步骤(2)的逻辑遍历,若可以找到下一个边缘点且两个连续边缘点所涉及的相同网格数<2则继续遍历,否则将VL从 中移除,此时若N(VR)<4,以VR为起点进行重新遍历。
(4)持续遍历到无法获取下一个边缘点为止,输出该矢量线。
从剩余的顶点中继续重复上述步骤,获取所有的边缘矢量线,这些矢量线组成了三维水质模型数据的外边缘和内边缘。
1.5 数字孪生网格体构建与推演仿真
1.5.1 水质孪生网格体构建与多模式呈现
静态网格体(Static Mesh)是数字孪生引擎中用于构建三维模型、地形场景等的重要数据类型,具备优异的渲染性能和效果,但需预先生成且无法动态更新。UE5中的程序化网格体(Procedural Mesh)克服了上述缺点,可以在不改变网格拓扑关系的前提下进行快速构建与更新。由于水质网格的拓扑关系在推演过程中维持不变,仅需更新网格顶点的高程、渲染颜色等属性,因此可利用程序化网格体数据结构,结合水质推演的动态信息,构建三维水质孪生网格体。
在数字孪生仿真引擎中,基于三维水质各分层的网格顶点和拓扑关系,以水质层为单元,分别创建水质分层和水质边缘的数字孪生网格体。对于水质分层,需要将污染物浓度依据分级色带转换为顶点渲染的颜色,在平面上平滑渲染为渐变颜色;对于水质边缘,需要将分层水深转换为网格顶点的高程,进行垂向分层封闭构造和法线计算,得到可自由控制的封闭边缘。
在完成三维水质数字孪生网格体的构建后,可进一步结合整体抬升、垂向夸张、边缘拉伸、纹理透明、分层显隐控制操作,满足三维水质数据的多样化展示需求,如图4所示。
1.5.2 水质指标演化过程时序推演仿真
根据三维水质网格对应的时序属性数据,对数字孪生场景中的水质分层网格进行重新渲染,在不改变数字孪生网格体拓扑结构的条件下快速更新显示效果,生成随时间变化的水质推演仿真结果。在推演时间发生变化的同时,同步更新数字孪生场景中环境天气效果、站点水位流量、工程设施运用等信息。
2 应用实例
2.1 研究区域概况
丹江口水库是南水北调中线工程水源地、国家一级水源保护区,水质安全保障尤为重要。丹江口水库水域面积广、深度大,污染来源分散且类型多,水库水质保护治理是一项复杂的系统工程,管理难度极大。2021年10月10日,丹江口水库首次蓄水至170 m,其后水库一直保持高水位运行。与原157 m正常蓄水位相比,水库淹没面积增加305 km2。环库存在大量消落带土地首次浸没。
2021年8—10月,丹江口水库入库、出库月平均流量较上年同期分别上升315.1%和301.5%,汉江陕西出境流量占入库流量近70%,2021年羊尾断面总磷浓度同比升高180%(2020年总磷浓度均值0.025 mg/L,2021年均值0.071 mg/L),带来了大量的污染物负荷,对库区水质产生不利影响。根据历史数据分析,丹江口水库总磷浓度10 a来呈逐渐升高态势,尤其是2021年丹江口水库总磷年均浓度达到0.023 mg/L,是近10 a来的最高值,已经达到Ⅱ类水质标准值上限(0.025 mg/L)的92%,存在较大的超标风险。
为了准确分析和模拟丹江口水库2021年秋汛期间总磷分布推演对陶岔供水安全的影响,需要开展模型计算与复演分析。
2.2 三维水动力水质模型计算
2021年秋汛期间,丹江口库区上游大量面源污染汇入水库,水库总磷浓度大幅升高。为保障陶岔供水安全,需要采用优化的调度方案,将丹江口水库下泄流量增加到10 878 m3/s,考虑在该方案下,水位最大值为170 m。
以库区上游水质监测站的总磷浓度实测值、水文站的水位和流量实测值、优化调度方案的水库的下泄流量等条件为输入参数,调用三维水动力水质模型,复演磷污染物入库后的动态演进过程;并以此为基础进行平面网格和垂向网格的划分,以12 h为推演输出时间间隔,计算展示8月21日至9月10日共计21 d的三维水质时空演变仿真结果。
在平面网格划分方向上,对丹江口库区采用正交四边形网格建模,如图5所示,共划分了17 261个网格,其中选用了2处有代表性的区域进行局部放大展示,局部区域1湖中间有小岛,局部区域2为河道转角处。
在垂向网格划分方向上,丹江口库区河底高程最小值105.3 m,最大值164.1 m,在计算过程中,考虑水位最大值为170 m。网格划分使用Sigma-Zed (SGZ)坐标系以减少垂直梯度误差,并采用如图6所示的划分方式,每个平面网格的垂向分层数为N(N≤5),最深的区域划分为5层,前N-1层每层高度为12.94 m。
经模型计算后,得出8月21日至9月10日共计21 d的水质推演数据,以二进制文件的方式存储,包括水质网格平面坐标文件、河底高程文件、时序水深文件、时序流场文件、时序总磷浓度文件。
2.3 数字孪生场景下的三维水质推演仿真
为了实现三维水质的推演仿真,基于Unreal Engine 5游戏引擎和Cesium For Unreal插件进行扩展开发,搭建了数字孪生仿真平台,在此之上运用高精度地理空间数据、水文气象数据、水质专题数据、防洪调度规则等建立丹江口区域的数据底板,并结合水质推演业务提供场景信息交互功能,开发形成具有交互推演能力的丹江口数字孪生系统,在系统的水质安全板块中实现了三维水质动态推演仿真。丹江口库区三维水质计算结果在系统中的集成如图7所示。
由于丹江口库区范围较大,图8展示了在数字孪生场景中采用背景消隐、水质网格体整体抬升、垂向夸张、
边缘拉伸等可视化方法提升三维水质的表达效果,并截取了不同推演时刻的三维水质总磷浓度状况,通过时序结果图呈现三维水质的演化过程。其中时间T分别选取1、3、6、10、15 d,三维水质总磷的浓度范围为0.01~0.13 mg/L。当总磷浓度为0.01 mg/L时使用蓝色进行渲染展示,表示水质良好;总磷浓度为0.13 mg/L时使用红色进行渲染展示,表示水质较差。从时序水质推演结果中可以看出,洪水入库后高浓度水体从汉库向坝前推进,在优化的调度方案下,水库下泄流量到达了10 878 m3/s后,高浓度水体基本不会进入丹江口库区,陶岔水质等级稳定在地表水水质Ⅱ类标准以上(总磷浓度≤0.025 mg/L),有效保障了陶岔供水安全。
2.4 推演仿真结果分析
相较于传统的水质可视化方法,本文方法具有以下优势:
(1)通过对三维水质网格的顶点复用计算、网格边缘的提取和分层构建,减少了50%的水质网格的顶点数量,减少了70%的边缘面片数量,降低了时序推演过程中的可视化模型更新复杂度,显著提升时序水质推演的帧率。
(2)以Unreal Engine 5游戏引擎+Cesium for Unreal插件为基础进行扩展开发,将三维水质数据按照真实的空间坐标转换为水质孪生网格体,可与高精度地形、实景三维模型等孪生场景要素无缝融合,并实现与水质推演区域测站、工程设施的交互查询与调度。
(3)将三维水质各分层信息构建为数字孪生网格体,提供水质网格抬升、垂向夸张、边缘拉伸、纹理透明、分层显隐等多种显示和控制模式,展示方法更加灵活、信息承载更加丰富,能直观、准确呈现分层水质指标,并渲染出连续、均匀过度的色彩效果。
3 结束语
本文提出了面向数字孪生场景的三维水质推演仿真方法,实现了三维水质推演过程的直观、准确、流畅、动态呈现,并以2021年丹江口水库秋汛总磷浓度分布推演为例,验证了方法的有效性。
下一步研究可继续深化和提升三维水质的推演仿真效果,包括:三维水质流场的动态推演仿真、主干流场趋势线的提取和可视化、典型水质断面的剖面提取和聚焦展示等方法。