《黄河数字流域模型.pdf》由会员分享,可在线阅读,更多相关《黄河数字流域模型.pdf(8页珍藏版)》请在得力文库 - 分享文档赚钱的网站上搜索。
1、第2卷 第7期 2007 年 7 月 492 黄河数字流域模型 王光谦,李铁键(清华大学水沙科学与水利水电工程国家重点实验室,北京 100084)摘 要:中国严重的水问题促使研究者在以数字化手段描述流域的同时,希望获得对流域中各种现象的模拟与预测能力。特别是在黄河流域,水资源短缺和水土流失是需要面对的两个主要问题。本文介绍了一个适用于黄河流域的基于物理机理的分布式的水沙过程连续模拟模型黄河数字流域模型。论文首先介绍数字流域的主要支撑技术,包括数字高程模型、流域分级与河网编码方法、模型参数提取和并行计算机制;其次介绍与黄河多沙粗沙区的水沙过程相对应的四个模拟模型,即坡面产流模型、坡面侵蚀产沙模型
2、、沟坡区重力侵蚀产沙模型和沟道不平衡输沙模型。各模型在数字流域平台上集成,形成流域尺度的水沙过程模拟功能。论文最后给出黄河数字流域模型在多沙粗沙区水沙过程模拟中的应用实例。关键词:数字流域模型;黄河;多沙粗沙区;流域泥沙;模拟 中图分类号:TV212-4 文献标识码:A 文章编号:16737180(2007)0704928 0 引言 在中国黄河中游的黄土高原,严重的土壤侵蚀与水土流失伴随暴雨-洪水过程频繁发生,一直以来都是造成中游土壤贫瘠、环境恶化,下游淤积严重、水患灾害频发的主要原因,从而使这一区域成为黄河水沙研究的重点。要对黄河流域进行全面模拟,就需要除水流之外,综合考虑流域内泥沙的产生、
3、输移和沉积的全过程,建立起基于物理机理的、能对全流域、全过程的水沙进行模拟与分析的系统。本文即为在清华大学的数字流域模型1,2平台上实现针对黄河的流域尺度的水沙过程模拟,即黄河数字流域模型的介绍。数字流域模型主题数据 后处理系统 模型层 地形源数据 遥感源数据 降雨量源数据 地理信息及数据提取系统 水保规划 减水减沙效益分析 灾害预警 应用层 数据层 下垫面源数据 坡面产流模型 坡面产沙模型 坡面模型:数据分析 数据可视化(GIS及 VR 技术等)数据挖掘 沟道汇流模型 沟道模型:重力侵蚀模型 河道动力学模型 图 1 数字流域模型的基本结构基金项目:国家自然科学创新研究群体科学基金(50221
4、903)通讯作者:E-mail: 中国科技论文在线 SCIENCEPAPER ONLINE 第2卷 第7期 2007 年 7 月 493 中国科技论文在线 SCIENCEPAPER ONLINE 数字流域模型定位于大范围、流域级的水与迁移物过程模拟,是一个具有多层空间分辨率的、模型参数易于获取的、能够实现并行计算的整体模型,其基本结构如图 1 所示。以数字高程模型、流域分级与河网编码方法、参数提取系统和并行计算机制等构成了其支撑平台,为模型层的实现和运行提供支持。为使作为核心的模型层符合水与迁移物运动的自然规律,按照自然过程的物理机理,将模拟模型分为坡面和沟道两类。其中坡面模型用于对坡面产流和
5、侵蚀过程进行模拟;沟道模型完成沟道汇流和污染物迁移的计算。图 2 黄土高原区的典型坡-沟系统3 图 3 黄土坡面的形态概化和模型概化示意图 而在黄土高原的主要水沙源区,即黄土丘陵沟壑区和高塬沟壑区,其地形支离破碎,植被稀少,土壤下渗能力低,水沙过程呈现出一些独有的特点。坡面的上端为地势平缓的峁顶或塬,此部分主要发生溅蚀。峁顶水流漫溢,到峁坡逐渐汇集,并冲刷形成细沟和浅沟。峁坡与沟道之间为沟坡区,坡面破碎,坡度极陡。峁坡上的浅沟往往通过跌坎与沟坡区的切沟相连,而切沟同时紧临沟道。沟道则位于沟坡的下部,是水流和泥沙的输运通道。黄土沟壑区的坡面-沟道系统(图 2)的单元概化如图3 所示。虽然在坡面-
6、沟道系统的不同部位发生着不同的侵蚀过程,但它是一个相互联系的整体,构成了侵蚀、输运、淤积等各种复杂的物理过程,并最终将泥沙带到流域出口。可概化出四个主要水沙运动过程,为:坡面产流、坡面侵蚀产沙、沟坡区重力侵蚀产沙和沟道系统的不平衡输沙过程。综上,可将与黄河数字流域模型相关的研究分为两类,一类是对流域信息的数字化描述和对模拟模型的支持技术,构成数字流域的支撑平台;另一类即为模拟模型,用以对流域内自然水沙过程进行模拟和预测。1 数字流域模型的支撑平台 构成数字流域支撑平台的主要技术包括数据库、GIS、RS、VR 和并行计算等内容,其中使数字流域模型不同于其他信息化领域的主要方面有:数字高程模型、流
7、域分级与河网编码、模型参数提取系统和并行计算机制,分述如下。1.1 数字高程模型 作为一种重要的地形数据形式,数字高程模型(Digital Elevation Modal,简称 DEM)是数字流域的基础,专业模型中所涉及的沟道密度,坡面坡度,沟道长度等都可以通过 DEM 提取出来。DEM 的出现使得具有物理机理的分布式水文模型成为可能,同时 DEM 也是产沙模型从经验性模型向物理机理性模型发展所必需的数据基础。但随着研究范围的扩大,大区域、高精度的 DEM 数据需要分块存储和读取,刘家宏等4进行了相关技术的研究,实现了从地形源数据获取、预处理到分块存储、读取拼接这一系列过程。DEM 的利用方式
8、有两种。一种将流域内的水分和迁移物运动方程直接在地形栅格上离散求解;另一种将流域划分为不同形式的子单元,假定各子单元内具有较为一致的地形和参数,应用积分形式的或简化的方程求解。数字流域模型采用了后一种方式,以坡面-沟(河)道作为基本单元(图 4),不同模型分别在坡面单元和由沟道组成的河网中进行计算。这样的结构与流域的自然过程具有一致性,便于按照自然过程建立具有物理机理的专业模型。第2卷 第7期 2007 年 7 月 494 图 4 坡面-沟道单元示意图 坡面和河网的提取一般采用 D85及其改进方法。首先利用DEM 计算得到流向,然后根据流向分析得到流域的河网、边界以及子流域的划分情况。数字流域
9、模型采用了基于此类方法的TOPAZ 模块6。1.2 流域分级与河网编码 专业模型的计算需要按汇流顺序在众多坡面-沟道单元上先后进行,并最终得到全流域的模拟结果。Strahler 级别能够表示汇流的先后顺序,但为了实现对流域由粗到细的分级描述和河网上高效的拓扑关系运算,数字流域模型实现了相应的河网编码方法。根据多数河流的树状形态特征,李铁键等7提出了基于二叉树的二元形式河网编码,能够表示河网上下游的连接关系和干支流的主次关系,并可直接索引。图 5 数字流域模型多沙粗沙区分区编码示意图 为适应划分子流域的习惯,对于覆盖范围广、地形精度高的大规模河网,采用了分级分区编码的方法。每一个分区(子流域)内
10、直接应用上述河网编码,而不同的分区具有不同的河流级别和位置。河流级别参照 Stralher 提出的地貌几何定量数学模型分级方法,将干流定义为第 1 级,支流的级别依次增长。同级支流按汇入干流的位置编号,最接近干流出口的编为第 0 号,向上游递增。将数字流域的河网编码方法应用于黄土高原多沙粗沙区的主要部分,划分为 4 级 37 个子流域,将各子流域不同级别的位置编号占两位数字按级别顺序连接,得到多沙粗沙区的流域分区如图所示,共完成了对 84618 个坡面-沟道单元的编码。1.3 流域单元的参数提取 数字流域模型和以往集总式流域模型的区别之一就是参数来源不同。数字流域模型以数字地形数据为基础,并综
11、合应用遥感技术,其大部分参数是随空间甚至时间分布的。各专业模型需要的坡面面积、坡面坡度、河段长度等几何信息可以通过地形数据计算获得,并将其对应到相应的单元存储。遥感图像是数字流域下垫面参数的主要来源,各专业模型所需要的植被覆盖、土地利用类型、流域土壤侵蚀系数、蒸发能力等都可以通过同时期的遥感图像解译或相应的专题图得到。为将图像中的信息对应到每个模拟单元,通过后者的中心点坐标或其多边形范围来获取对应专题图上的像素值,经过转化而得到相应参数。通过数字地形和遥感图像获得的参数储存在数据库中,以提供给专业模型使用,其中每一条记录对应于一个坡面-沟道单元。1.4 并行计算机制 成熟的集群硬件技术使各专业
12、领域的并行算法研 究 和 应 用 层 出 不 穷。MPI(Message Passing Interface,消息传递接口)8是目前世界上最重要,也最流行的并行编程标准。MPI 通过独立于计算机语言的函数库实现,它具有移植性好、功能强大、效率高等多种优点,而且有多种免费的实现版本,支持不同的软硬件平台。其 MPI-2 版本提供了 287 个用于通讯与管理等不同功能的函数,并在 C/C+和Fortran77/90 语言下有稳定高效的实现。在数字流域平台上进行的水沙过程等专业子模型计算带来了巨大的计算量,如果采用串行算法,将会耗费大量时间,且不能充分发挥数据库的效率。数字流域模型的单元间具有显著的
13、数据弱相关的特点,非常符合并行计算的条件。并行化的方式是在大规模河网上的空间域分解。为处理水流演进和物黄河数字流域模型 第2卷 第7期 2007 年 7 月 495 中国科技论文在线 SCIENCEPAPER ONLINE 质输送的上下游依赖关系,采用了动态任务分配的形式,并为提高效率而提出了一种并行粒度的动态控制方法。模型的并行计算也因而采用主从模式的组织形式,由专门的控制节点完成任务调度9。2 黄河数字流域的模拟模型 数字流域模型的支撑平台是基础,模拟模型是核心。黄河数字流域模型建立了以下模拟模型,用于完成流域水沙过程计算,具体内容可参考文中给出的相关文献。2.1 坡面产流模型 图 6 坡
14、面的剖面结构图 WdfWufPEcan Qs Qgu Qgd Eu河网ScanWuWdPnqzuqzdWusSc0Wds植被层表层土壤下层土壤 图 7 产流模型流程示意图 黄土高原土层包气带可深达几十米,因此,本文忽略潜水对黄土坡面产流的影响,采用如图 6 所示的坡面结构,将坡面土体划分为表层土和下层土两部分。产流模型被设计为一个以描述超渗产流为主的机理型模型,在物理图景概化和产流机理符合黄土高原坡面产流自然过程的前提下,采用尽可能简单的形式模拟不同子过程。模型主要计算植被截流、地表超渗产流和表层土快速壤中流等子过程,并兼容可能出现的下层土壤中流和表层土短时超蓄产流。为实现多场次降水的连续模拟
15、,计算蒸散发过程和土壤内的水分迁移。计算流程如图 7 所示。冠层蓄水、表层土壤水和下层土壤水的水量守衡关系为:canncanEPPtS=;()guuzdzuuQEqqAtW=;gdzddQqAtW=(1)式中 Scan为冠层蓄水量(单位 m),t 为时间(单位s),P 为降雨强度(单位 m/s),Pn为穿过植被冠层的净雨强(单位 m/s),Ecan为冠层积水的蒸发强度(单位m/s);Wu为表层土蓄水量(单位 m3),A 为坡面的投影面积(单位 m2),qzu为地表下渗速率(单位 m/s),qzd为表层土与下层土间水分垂向运动速率(单位 m/s,向下为正),Eu为表层土壤蒸发强度(单位 m/s)
16、,Qgu为表层土壤出流量(单位 m3/s);Wd为下层土蓄水量(单位 m3),Qgd为下层土壤出流量(单位 m3/s)。其中植被冠层截流计算采用 MIKE SHE WM 中Kristensen-Jensen 模型10,11的形式。表层土壤蒸发量采用傅抱璞提出的公式12计算,下层土参与蒸发的过程通过其中水分向上层土的调整间接体现。地表超渗产流和两层土间的水分调整基于土壤水动力学和土水势理论计算。当土层含水量超出田间持水量时,计算壤中水出流。2.2 坡面侵蚀产沙模型 产沙模型将溅蚀、片蚀和沟蚀过程统一为侵蚀量沿坡面增加的概化过程13。将侵蚀量与水流条件挂沟,通过对坡面水沙运动的连续条件进行推导,得
17、到坡面单位时间的侵蚀量为:32332555105542sEmDA Q nJL=(2)式中 E 为坡面总侵蚀量(单位 kg/s);为泥沙与水流运动速度的比值;A 为坡面面积(单位 m2);m为单位坡长下的侵蚀层数增量(单位 m-1);Q 为坡面流量(单位 m3/s);n 为曼宁系数;J 为坡面坡度;L 为坡面长度(单位 m)。m 值与水流侵蚀能力和土壤抗冲性有关,可写作 33755100mmkkqnJ=(3)第2卷 第7期 2007 年 7 月 496 式中 k、为与土地类型和土壤抗冲性有关的系数,在计算前率定,且1。将式(3)代入式(2)得到坡面侵蚀产沙量的最终表达式为:23332735155
18、551010542sQEAkg qnL DJA+=(4)2.3 重力侵蚀模型 重力侵蚀建立在沟坡区,以将发生滑坡或崩塌的土体为对象,考虑水流的诱发作用,对土体进行力学分析14。如图 8,其受力状况包括:1)土体重力,并考虑由于水流入渗产生的重力增量;2)滑裂面上的正应力与切应力,并考虑由于水流入渗而引起的粘聚力降低;3)滑坡体顶端沿黄土垂直节理的初始裂缝产生的抗滑力减少和正向水压力;4)沟道水流切割导致的坡脚后退或湿陷。沟道水流径流坡脚淘刷坡脚湿陷坡脚淘刷坡脚湿陷水流入渗G增大增大减小减小雨降T 图 8 重力侵蚀土体受力示意图 其中,渗流过程采用产流模型的计算结果。土体强度采用土力学方法计算,
19、并随含水量变化。坡脚冲刷采用 Osman 的方法15,土体单位时间侧向冲刷距离为:()1.3clcsCeB=(5)式中B为土体单位时间受水流侧向冲刷而后退的距离(单位 m);Cl为与土体理化特性有关的系数,根据试验资料,可取43.64 10lC=;为水流切应力(单位 Pa);c为土体起动切应力(单位Pa)。将土体受力归为下滑力与抗滑力两类。可分别表示为:sincosDtFWT=+(6)tanNcLFR+=(7)式中 FD为土体下滑力;Wt为土体重力;为滑裂面角度;T 为土体顶端裂隙水压力。FR为土体抗滑力;L 为滑裂面长度;N 为滑裂面上的正应力。土体这些受力均可以通过其几何形态和土的力学性质
20、计算得到,并在时间上随含水量而变化。通过对土体的下滑力和抗滑力进行分析,可以得到其稳定安全系数。为了符合重力侵蚀发生的随机性,对安全系数进行模糊概率分析,得到其失稳隶属度,并在模型运行时对随机事件的发生与否进行判断,最终完成对重力侵蚀的模拟。2.4 沟道系统水沙运动模拟 沟道系统中水沙运动的模拟是依上下游顺序逐河段完成的,由于泥沙源区的沟道不具备实测断面资料,仅有河段长度、河底平均坡降等基本参数,因此水流模拟采用以扩散波方法为基础的形式,并假定沟道横断面形状为“V”型,以获得水深、流速等水力参数用于输沙计算。扩散波方程系数采用四点平均格式计算:4=iCC,()CQfB,=,()4/=iiCQC
21、Q;4,3,2,1=i (8)式中 C 为波速(单位 m/s),B 为水面宽(单位 m),表示计算点的值,其中C用于计算槽蓄系数K,B和CQ用于计算比重因子。考虑水面附加坡降,取水力坡降 S 为:xhSS+=0 (9)式中 S0为河床坡降,为系数,取值范围为:0.7 1.0,0,=涨水时落水时 (10)黄河数字流域模型 第2卷 第7期 2007 年 7 月 497 中国科技论文在线 SCIENCEPAPER ONLINE 将xh转化为已知量Q的函数:xQBCxhh=1 (11)式中hQBCh=1,为水深形式扩散波方程的波速系数(单位 m/s)。为控制负反应,需使空间步长:+BSCQtCBSCQ
22、tCx,(12)因而在计算中时间步长先确定的情况下设置每个河段的分段数:=tCLIntN,1N (13)式中L为河段的长度(m),()Int为取整运算。输沙计算采用悬移质不平衡输沙的模式。对于输沙能力,可采用张红武公式或费祥俊提出的沟道高含沙公式16两种方法计算。3 在多沙粗沙区水沙过程模拟中的应用 模拟以 1977 年为例。1977 年虽然属于中等水量年份(年降雨量 513mm),但是局部地区出现了罕见的特大暴雨,特别是孤山川、延水和无定河,而该区域一直以来又是侵蚀产沙最为严重的区域。因此,该年不仅输沙量大,而且输沙强度高,具有很强的代表性。几个主要支流把口站的计算含沙量过程与实测过程的对比
23、如图 9图 12 所示,计算与实测输沙量统计对比见表 1。更进一步,得到多沙粗沙区出口站龙门的模拟水沙过程如图 13 和图 14 所示。其中头道拐入流采用的是实测数据。0200400600800100012007月2日7月12日7月22日8月1日8月11日8月21日8月31日观测日期含沙量(kg/m3)实测含沙量计算含沙量 图 9 皇甫川皇甫站含沙量过程计算与实测对比 01002003004005006007007月2日7月12日7月22日8月1日8月11日8月21日8月31日观测日期含沙量(kg/m3)实测含沙量计算含沙量 图 10 孤山川高石崖站含沙量过程计算与实测对比 010020030
24、04005006007008007月2日7月12日7月22日8月1日8月11日8月21日8月31日观测日期含沙量(kg/m3)实测含沙量计算含沙量 图 11 窟野河温家川站含沙量过程计算与实测对比 01002003004005006007008009007月2日7月12日7月22日8月1日8月11日8月21日8月31日观测日期含沙量(kg/m3)实测含沙量计算含沙量 图 12 佳芦河申家湾站含沙量过程计算与实测对比 01000200030004000500060007000800090007月1日7月11日7月21日7月31日8月10日8月20日8月30日观测日期流量(m3/s)计算流量实测流
25、量 图 13 龙门站计算流量过程与实测值对比 第2卷 第7期 2007 年 7 月 498 表 1 1977 年黄土高原多沙粗沙区主要支流年输沙量 支流名称 实测输沙量(108t)计算输沙量(108t)误差()皇甫川 0.26 0.33 26.92 孤山川 0.839 1.21 44.22 窟野河 1.38 1.52 10.14 秃尾河 0.211 0.18-14.69 佳芦河 0.121 0.17 40.5 三川河 0.465 0.69 48.38 无定河 2.69 3.42 27.14 清涧河 1.17 0.95-18.80 05001000150020002500300035004000
26、450050007月1日7月11日7月21日7月31日8月10日8月20日8月30日观测日期输沙率(t/s)实测输沙率计算输沙率 图 14 龙门站计算输沙率过程与实测值对比 图14龙门站计算输沙率过程与实测值对比从年输沙量来看,多沙粗沙区从北部的皇甫川到南部的清涧河之间的 8 条河流,计算得到的年输沙量与实测值相比,数量级是一致的。从输沙过程来看,计算日均含沙量过程与实测过程基本趋势是一致的,能够反映多沙粗沙区产沙频繁,含沙量变化幅度较大的特点。总的来说,这些结果基本反映了汛期多沙粗沙区对黄河中下游的水沙输入情况,并提供了龙门站的模拟水沙过程,能够完成对多沙粗沙区水沙过程的预报。4 结论与展望
27、 黄河数字流域模型借助先进的信息技术手段,在支撑平台上整合坡面产流产沙、沟坡区重力侵蚀、沟道水沙运动等多个与流域水沙过程相对应的模拟模型,能够在流域尺度上进行分布式的水沙过程全面模拟。应用实例表明模型具有可以接受的精度。数字流域作为一种平台应能够承担更多内容和形式的应用。今后的工作一方面在于增强和改进各模拟模型的机理与计算精度,另一方面要针对不同的区域和应用领域,建立相应的应用系统,为流域的可持续发展提供相适应的支撑。参考文献 1 王光谦,刘家宏.数字流域模型M.北京:科学出版社,2006.2 王光谦,刘家宏,李铁键.黄河数字流域模型原理J.应用基础与工程科学学报,2005,13(1):18.
28、3 张天曾.黄土高原论纲M.北京:中国环境科学出版社,1993.4 刘家宏,王光谦,王开.大流域 DEM 数据存取关键技术研究 J.清 华 大 学 学 报(自 然 科 学 版),2004,44(12):16461649.5 OCallaghan J F,Mark D M.The extraction of drainage networks from digital elevation data J.Computer vision,Graphics and image processing,1984,28:323344.6 Garbrecht J,Martz L W.TOPAZ overvie
29、w Z.USDA-ARS,Grazinglands Research Laboratory,1999.7 李铁键,王光谦,刘家宏.数字流域河网编码方法J.水科学进展.2006,17(5):658664.8 Message Passing Interface ForumEB/OL,http:/www.mpi-forum.org/index.html.2006.9 李铁键,刘家宏,和杨,等.集群计算在数字流域模型中的应用J.水科学进展.2006,17(6):869-874.10 Vzquez R F,Feyen J.Effect of potential evapotranspiration es
30、timates on effective parameters and performance of the MIKE SHE-code applied to a medium-size catchment J.Journal of Hydrology.2003,270(3-4):309327.11 Kristensen K J and Jensen S E.A model for estimating actual evapotranspiration from potential evapotranspiration J.Nordic Hydrology.1975,(6):7088.12
31、傅抱璞.土壤蒸发的计算J.气象学报.1981,39(2):226236.13 王光谦,薛海,刘家宏.坡面产沙理论模型J.应用基础与工程科学学报,2005(增刊):1-7.14 王光谦,薛海,李铁键.黄土高原沟坡重力侵蚀的理论模型黄河数字流域模型 第2卷 第7期 2007 年 7 月 499 中国科技论文在线 SCIENCEPAPER ONLINE J.应用基础与工程科学学报,2005,13(4):335-344.15 Osman A.M,Thorne C.R.Riverbank stability analysis I:theory J.Journal of Hydraulic Enginee
32、ring,1988,114(2):134150.16 费祥俊,邵学军.泥沙源区沟道输沙能力的计算方法J.泥沙研究,2004,1:18.Digital Yellow River Model Wang Guangqian,Li Tiejian(State Key Laboratory of Hydroscience and Engineering,Tsinghua University,Beijing 100084)Abstract:Soil erosion is one of the key concerns in land use management for the Loess Platea
33、u of the Yellow River,where serious soil loss is the root cause of environmental and ecological degradation of the basin.In this paper,a physically-based,distributed-parameter,and continuous erosion prediction model at the river basin scale was developed with the aim of assisting in developing bette
34、r land use management strategies.The framework,the major supporting techniques,and the typical erosion processes are described.The physical processes of sediment yield and transport in the Loess Plateau are divided into four sub-processes,including the runoff yield on hillslopes,hillslopes soil eros
35、ion,gravitational erosion in gullies,and hyperconcentrated flow routing in channels.For each sub-process,a physically-based simulation model was developed and embedded into the whole model system.The model system was applied to simulate the runoff and sediment yield in several typical years in the coarse sediment source area of the Loess Plateau,and the simulated results were in reasonably good agreement with the measured values.Key words:Digital Watershed Model;the Yellow River;soil erosion;sediment yield;simulation
限制150内