（水文学及水资源专业论文）矩阵标识法在感潮河网非恒定流计算中的应用.pdf

矩阵标识法在感潮河网非恒定流计算中的应川 水文学及水资源 廖梓瑾 指导教师：刘美南副教授 摘要 三角洲平原河网地区通常是城市发达、人口众多的地区，同时又是水域密布、 地势较低，易于发生洪涝灾害的地区。近年来，由于自然和人为多方面的因素影 响，引起了平原河网水文条件的复杂多变以及水环境质量的日益恶化，成为影响 经济可持续发展的一个主要问题。为了模拟河网非恒定流运动规律，本文在前人 研究的基础上，建立了三角洲河网水动力计算模型。应用三级联解法的思路，采 用改进的迭代法一一矩阵标识法求解连接矩阵。该法克服了以往方法依赖汉点编 码优化的缺陷从而避免了河网计算中微小的河道改动导致重复优化编码的不便： 具有编程简单、计算效率高，可移植及通用性强的特点。珠江三角洲河道众多， 纵横交错，给河网的概化带来困难，在模拟过程前期需对概化河网进行修改，因 此将矩阵标识法应用于珠江三角洲河网的水动力模拟，可极大地提高工作效率。 取典型枯水期1 9 9 1 年1 2 月1 5 日的实测水文资料对5 个水位站点进行参数率定及 模型验证，以率定好的模型对1 2 月1 6 目至1 7 日的河网水流进行计算，取得了较 好的模拟效果。对河网水动力场的准确模拟是合理进行平原河网的灌溉、排涝、 防洪、规划和污染控制等工作的基础，对指导具体生产实践具有重要意义。 关键词：河网非恒定流，三级联解，矩阵标识法，珠江三角洲 s t u d y o fm a t r i xm a r km e t h o di nt i d a lr i v e rn e t w o r k c o m p u t a t i o n h y d r o l o g ya n dw a t e ri _ e s o u r c e s l i a o z i - j i n s u p e r v i s o r ：p r o f e s s o rl i u m e i n a i l a b s t r a c t t h ed e l t ar i v e rn e t w o r ka r e a so f t e ni n c l u d em a n y d e v e l o p e d c i t i e sa n d p e o p l e ，a n d t h e a r e a sh a v es om a n yw a t e ra r e a st h a ti t s e a s yt ol e a dt o f l o o do rw a t e r l o g f o rt h e s e y e a r s ，t h eh y d r o l o g i cc o n d i t i o nh a sb e e nc o m p l i c a t e da n dt h ew a t e rq u a l i t yh a v eb e e n d e t e r i o r a t e db e c a u s eo f m a n y f a c t o r sa b o u tb o t hn a t u r ea n d p e o p l e ，w h i c h h a sb e c o m ea m a j o rp r o b l e me f f e c t i n gt h ee c o n o m i cc o n t i n u a n c ed e v e l o p m e n t f o rs i m u l a t i n gt h e w a t e rr u l ei nr i v e rn e t w o r k ，t h i sa r t i c l ee s t a b l i s h e st h ed e l t ar i v e rn e t w o r kh y d r o l o g i c m o d e l b a s i n g o nt h e f o r m e r s s t u d y i n g t h e m o d e lu s e st h em e t h o dc a l l e d c h a n n e l - j u n c t i o n c h a n n e l a n dm a t r i xm a r k m e t h o d ，w h i c h c a no v e r c o m et h e d i s a d v a n t a g e so fs a v i n gc o m p u t e rs t o r a g ea n dr a i s i n gc a l c u l a t i o ne f f i c i e n c y b e s i d e s ，i t i su n r e s t r i c t e di nn o d ec o d i n g ，a n di th a s g o o de x p a n s i b i l i t y , t r a n s p o r t a b i l i t y a n d u n i v e r s a l i t y t h ep e a r lr i v e r d e l t ah a sl o t so fr i v e r sa n di t sd i f f i c u l tt os i m p l i f yt h er i v e r n e t w o r k s ou s i n gt h em a t r i xm a r km e t h o dt ot h ep e a dr i v e rd e l t ac a ni m p r o v e w o r k i n ge f f i c i e n c y c h o o s i n gt h et y p i c a ll o w w a t e rd a t a ，d e c 1 5 ，1 9 9 1 ，t h em o d e lg e t s p a r a m e t e r sf r o mf i v ew a t e r l e v e ls t a t i o n s t h e nt h em o d e lc o m p u t e st h eh y d r o l o g yf r o m d e c 1 6t o 1 7 ，a n dg e t sg o o de f f e c t t h es i m u l a t i o no fw a t e rf l o wi s t h eb a s eo f i r r i g a t i o n ，e j e c t i n gw a t e r l o g g i n g , p r e v e n t i n gf l o o d ，p r o g r a m m i n g a n d c o n t r o l l i n g p o l l u t i o n ，e t c k e y w o r d s ：d e l t ar i v e rn e t w o r k ，c h a n n e l - j u n c t i o n c h a n n e l ，m a t r i xm a r k m e t h o d ，p e a r l r i v e rd e l t a 矩阵标识法枉感潮河网非恒定流计算中的应j 】 第1 章绪论 1 1 选题背景和研究意义 在平原河网地区，存在大量分汉河道，且由于多存在于河流下游入海口，受 潮汐影响显著，多为非恒定流。 三角洲平原河网地区通常是城市发达、人口众多的地区，同时又是水域密布、 地势较低，易于发生洪涝灾害的地区。近年来，一方面由于全球性气候异常，洪 水、暴雨和强潮的综合影响引起了平原河网水文条件的复杂多变；另一方面，随 着城市人口的增多、经济的发展，各种工业废水和生活污水大量排入河网，水环 境质量日益恶化，成为影响经济可持续发展的一个主要问题。所有这些对平原河 网的灌溉、排涝、防洪、规划和污染控制等提m 了一系列新的研究课题。由于三 角洲河网区域范围很大，水流条件复杂多变，大多数情况下只能采用数值方法进 干亍模拟，其核心问题是河网数学模型的建立及求解，其中水动力模型又是其它模 型，如水质模型、水环境容量等模型的基础【l “。 珠江三角洲位于广东省陆域南部，是西江、北江和东江下游的冲积平原。范 围包括西江、北江汇合口思贤活以下的西、北江三角洲和东江石龙以下的东江三 角洲，以及入注三角洲诸河。珠江三角洲是广东省最大的平原区。三角洲内星罗 棋布的丘陵、台地、残丘使水道分八口门出海，所谓三江汇合、八口分流就成为 珠江三角洲的主要地貌特征。流域而积2 6 8 2 0 k m 2 ，占珠江全流域面积的5 9 1 【1 1 l 。 珠江三角洲地理位置和经济环境条件优越，是我幽的发达地区，并持续保持着较 高的经济增长速度，同时也是世界关注的几大经济发展区之一。但近1 0 年来，由 于多种因素的共同影响，三角洲洪涝灾害及威胁日趋加大，与区域经济发展态势 极不协调。据统计资料分析1 1 2 , t 3 l ，1 9 9 4 年6 月洪水期间，中山市涝渍水位1 7 m ， 直接经济损失7 5 3 亿元。1 9 9 8 年，由于洪水和台风暴潮，三角洲的洪涝灾害直接 矩阵标识法在感潮河喇非恒定流计算中的应用 经济损失就达4 2 5 亿元。同时，三角洲也常常是废污水及污染物的主要排放区域， 以珠江三角洲为例，其废污水排放量约占广东省的7 5 ，而河网区内废污水排放 量约占珠江三角洲的7 5 5 ，c o d 排放量约占6 8 1 ，氨氮排放量约占6 7 8 1 ”i 。 珠江三角洲河网区受潮汐影响显著，水系构成和水沙动力条件都十分复杂， 潮水和径流之间的相互关系年内变化复杂，口门区潮流多为往复流。潮汐河网区 水流受径流和潮流的共同作用，其水流动力学特征既不同于内陆顺流河，又有别 于海洋，是潮汐河口中最为复杂而且典型的一类水体。强烈的人类活动，复杂的 网河地形地貌，江海相接。多口门分流，水道纵横交错，会潮、会流点众多，构 成了潮汐河网水文水动力学上流态多变、水情复杂的特点：因此，为了更好地进 行灌溉、排涝、防洪、规划和污染控制等工作，需要建立通用性强、方便实用的 河网水动力计算模型来模拟水流运动规律。 1 2 浅水动力学计算研究现状 在流体力学中浅水流动是对实际流动的一种简化和概化数学模型。1 8 7 1 年法 国人圣维南f a j c b d e s a i n t v e n a n t ) 建立了后来以他名字命名的一维非恒定水流方 程组【1 5 】，其中包含了浅水流动的三个基本假设，即：( 1 ) 流速沿座标轴方向并在 横断面内均匀分布( 底摩阻的影响局限在边界层内) ；( 2 ) 水压力沿水深服从静压 分布；( 3 ) 无垂直速度与加速度，它们或是由流动的非一维性所引起( 如断面形状 及弯道产生的二次环流) ，或是由密度不均匀产生的浮力所引起( 如水中含沙、盐、 热、气等导致的垂向交换) 。虽然严格意义的浅水流动并不存在，它接近于但不等 同于某些实际流动。除深海外大部分地表水在一定条件下都可近似作为浅水水体 来处理。传统水力学的主要对象亦是用实验方法研究浅水流动。后人把圣维南方 程组( 可称为一维浅水方程) 推广为二维浅水方程组，只要把上列第一条浅水假 设改述为水平流速沿垂线均匀分布即可。 谭维炎在文献 1 5 q u n - 维n o - - 维浅水流的数值模拟进展历程和当代前沿问 题作了回顾。 2 矩阵标识法在感潮河网非恒定流计算中的应用 1 2 1 一维渐变浅水流的数值模拟 自5 0 年代首次应用计算机和差分法( f d m ) 模拟河道水流，至7 0 年代中出版 y e v j e v i c h 等主编、由当时一批著名专家分工执笔的明渠非恒定流三卷本名著， 标志了一维渐变浅水明渠流数值模拟的初步成熟。该书从极其大量文献中( 书目 构成该书第3 卷) 选出几种差分格式加以推荐。随着格式的不断发展成熟，其中 l a x w e n d r o f f 显格式已为计算气动力学1 9 6 9 年提出的m a c c o r m a c k 格式所代替( 二 者均属二阶中心格式，为计算方便而修改前者成为后一格式) ，蛙跳及扩散格式等 极少使用，p r e i s s m a n n 四点隐格式一直保持其强劲生命力，特征格式亦不断有所完 善。近年来，美、法、荷、丹、日等国家的著名国际工程顾问公司所使用的、及 投放国际应用软件市场的一维水流模拟软件，大多使用二阶中心f d m ( 包括箱格 式) 。其原因是：f d m 建立在经典的数学逼近理论的基础上，简单且易为接受， 处理效率较高。 长期以来一维浅水流算法不能用于模拟急变和间断流。自7 0 年代中起，计算 气动力学对激波模拟取得了突破性进展。在对可压流欧拉方程组进行深入研究的 基础上，先后提出了几类深刻的思路和巧妙的技术。其中绝大多数的目标都是试 图建立守恒逆风格式( 因为以特征分析为基础故称为特征基逆风) ，因为只有守恒 型格式才能计算间断( 当然，此时不守恒误差不存在) ，同时只有逆风格式才含足 够的数值粘性以消除虚假振荡( 中心格式粘性太小) 。一类是从代数特征分析出发 的通量向量分裂( f v s ) 格式。二类是使用黎曼解的格式，包括通量差分裂( f d s ) 格 式，及按特征分解的o s h e r 格式等。三类是在无论使用中心或逆风格式时，对其中 所含的反扩散项或粘性项加以限制，以达到既无虚假振荡又保持较高精度的目的， 其中包括7 0 年代初提出的通量校正输运( f o r ) 格式，及8 0 年代初提出的全变差消 减( t v d ) 格式和8 0 年代末的基本无振荡( e n o ) 格式等。此外，还有许多格式小类， 包括以m o r e t t i 为首的意大利学派的 格式( 传统特征法的发展) ，和法国学派的中 心格式加上隐式校正等。 3 卵眸标识法在感潮河网非恒定流计算中的应用 1 ，2 2 二维浅水流的数值模拟 1 9 6 7 年美国l e e n d e r t s e 首次应用交替方向隐( a d i ) 差分格式模拟二维潮汐水 流，并很快得到推广。后来国际上又应用苏联y a n e n k o 等人的分裂法，从按空间 座标分裂到按物理机制( 即方程中物理意义不同的项) 分裂控制方程组。我国在8 0 年代也经历了类似的发展过程。这些算法的共同点都是基于矩形座标网格上的 f d m 。明显的缺点是，受网格密度的限制，对平面水体的周边形状逼近较差，且 对边界条件的实现困难或不准确，尤其是边界线的转折点。 7 0 年代起计算流体力学界开始将有限元法( f e m ) 应用于平面流动模拟。其原 理是分单元对解逼近，使微分方程空间积分的加权残差极小化。由此建立f e m 方 程给出数值解。f e m 的最大优点是所用的不规则的无结构网格，易于较准确地逼 近浅水体周边地形和水下地形。但f e m 用于水流计算问题时，对非恒定流，每一 时间步要求解一个大型线性方程组，耗机时多。同时f e m 在数学上适于求解椭圆 型方程的边值问题。不适于求解以对流为主的输运问题，即使加以改进如逆风、 最小二乘等方法，结果也未能令人满意。再次，f e m 的性能类似中心差分格式， 缺乏足够耗散，不适于计算间断，要加上人工粘性。因此，f e m 迄今在流动计算 中尚未得到广泛应用。 另一种计算方法为林秉南先生于5 0 年代初提出的特征法m o c ，它与f d m 的 主要不同在于利用沿特征成立的特征方程，而不是利用普通空间坐标中的原始方 程。特征方程反映了双曲问题中信息沿特征传播的性质，因而算法符合水流的物 理机制；是一种更合理的逆风格式，具有优良的精度。不足之处是。特征方程常 为非散度形式，即非守恒形式，用差分法离散特征方程时会带来守恒误差a 当水 流沿程变化较大时，如存在底坡，非齐次项的计算较繁，且可能带来可观的误差。 m o c 也不能直接计算间断解。 用有限体积法f v m 模拟二维浅水流在9 0 年代已呈明显优势。和f e m 一样将 计算域划分成若干规则或不规则形状的单元或控制体。在计算出通过每个控制体 边界沿法向输入( 出) 的流量和动量通量后，对每个控制体分别进行水量和动量 平衡计算，便得到计算时段末各控制体平均水深和流速。因此，f v m 是对于推导 原始微分方程所用控制体途径的回归，与f d m 和f e m 的数值逼近相比其物理意 矩阵标识法在感潮河网非恒定流计算中的应用 义更直接明晰。对由一个或多个控制体组成的任意区域，以至整个计算域，都严 格满足物理守恒律，不存在守恒误差，并且能正确计算问断。设计f v m 格式的关 键在于如何计算跨控制体界面的通量。f v m 能像f e m 一样适用于任意的不规则 网格，且着眼于控制体上的逼近，具有守恒性；又能像m o c 格式一样具有以特征 为基础的逆风性。并且，在具有上述优良性能的同时，处理效率与f d m 相近，而 远高于f e m 。可以说，f v m 体现了f e m 的几何灵活性、m o c 的精度和f d m 的 效率和守恒性。实现它的关键在于如何自动生成适用的网格，如何选择f v m 的形 式，如何估计通量，以及如何提高精度等。 1 3 河网水动力计算研究进展 河网非恒定流水流不同于单一河流的特点在于河网水流的错综复杂性，以及 由此带来的方程组离散和求解上的困难，是多年来河网问题研究的一大难点。虽 然f v m 在二维浅水流模拟中呈现优势 1 6 - 2 1 】，但在河网的水动力模拟上，现在仍多 采用f d m 离散一维圣维南方程组进行模拟。 1 ，3 1 河网水动力计算模型 河网水动力模型大体可以分为河道一汉点一河道模型、单元划分模型、混合 模型以及人工神经网络模型4 类1 2 2 1 。 河道一汊点一河道模型的基本思想是：将河网中的每一河道视为单一河道， 其控制方程均为一维圣维南方程组；河道汇合处称为汊点，每个汉点处均应满足 水流连续性方程和能量守恒方程。求解由边界条件、圣维南方程组和汉点衔接方 程组成的闭合方程组，即可得到各河段内部断面的未知水力要素。河道一汉点一 河道模型是目前最为流行的河网模型。 单元划分模型的基本思想【2 3 】是：将水力特性相似、水位变化不大的某一片水 体概化为一个单元，单元间流量交换的是连接河道，其本身无调蓄作用。假定单 元问存在两种连接方式：在无水： 建筑物或障碍物( 不存在局部水头损失) 的情 s 矩阵标识法在感潮河删非恒定流计算中的应用 况下认为是河型连接；堰型连接存在局部水头损失，又可分为自由出流和淹没出 流两种形式。取单元几何中心的水位为单元代表水位，给出水位与水面面积关系。 将计算河网分解为一定数量的单元，再进行分组，然后确定各单元问的连接类型。 对每个单元给出微分形式的质量平衡方程，经有限差分法离散后得到以单元水位 为基本未知量的方程组，进而求解各单元的代表水位和单元间流量。 姚琪等人认为，河道一汉点一河道模型和单元划分模型都不能很好地适应运河 水网的特性，前者失之过繁，后者失之过简【川。混合模型是将河道一汉点一河道模 型和单元划分模型中与平原河网特性相适应的优点综合起来，并避免其不相适应的 缺点而构成的新模型。其建模基本思想是：将平原河网的水域区分为骨干河道和成 片水域两类，对骨干河道采用河道一汉点一河道模型：对成片水域采用单元划分的 方法将其划分为单元，再引入当量河宽的概念，把成片水域的调蓄作用概化为骨干 河道的滩地，将其纳入河道一汉点一河道模型一并计算。 人工神经网络模型是采用计算机来模拟生物神经网络的结构和功能，由大量神 经元构成的并行分布式系统。人工神经网络与平原河网在结构上有许多相似之处， 两者都是由各个内部单元通过并联或串联形成一个相互制约的整体网络结构，通过 调整系统内部各个“神经元”之间的相互作用可以达到系统输入变量和输出变量之 间的最优化或平衡，因此，人t 神经网络理论可用于复杂河网水动力模型的数值模 拟【删。但天然河网十分复杂，河流湖泊众多，很难直接应用神经网络模拟河网的各 种参数，必须首先根据实际工程问题研究的需要、河道湖泊之间的相互关系以及数 值计算的需要，对河网进行概化。李荣等人将人工神经网络的原理应用于荆江和洞 庭湖灌溉系统的水量预报i ”1 。 河道一汉点一河道模型原则上可以求解任何河网的水力参数，也是目前应用最 广的河网模型。 单元划分模型使河网水系得到明显的概化，这对于湖泊、池塘、水库众多的平 原复杂河网，可使模型的数值计算大为简化。其缺点是：水位相近是单元划分的主 要依据，忽略了主要河道、支流河道、湖泊其他水力特性的差别；忽略了圣维南方 程组的惯性项，仅适用于河道流速时空变化不大的情况，对于汛期洪水涨落比较急 剧和沿海感潮河段，该模型不适用。 混合模型综合了河道一汉点一河道模型和单元划分模型的优点，对平原复杂河 网进行合理的概化，进而归入河道一汉点一河道模型一并计算，可以大大简化河网 矩阵标识法在感潮河网非恒定漉计算中的应用 水流的计算。不足之处是河网前期概化工作量大，成片水域和骨干河道的划分属经 验处理，不易掌握。 人工神经网络模型和水动力模型的有机结合避免了水动力模型计算量大、计算 速度难以满足实时预报的要求等问题；同时可以利用水动力模型为神经网络模型提 供学习样本，弥补了在重要河段和区域由于资料缺乏而使神经网络模型应用受到局 限的缺陷。神经网络模型是“黑箱”模型，反映网络内部节点之间关系参数的物理 意义模糊，就河网水动力模拟而言，它实际上反映的是河网各种因素的综合影响， 这就给模型的验证带来困难，同时要求神经网络模型要有足够长的学习样本，否则 将影响计算结果的精度。文献 2 6 1 仅就河网水量的数值模拟建立了神经网络模型， 未涉及到河网水流运动的另外一个重要参量一一水位。 1 3 2 对河道一汉点一河道模型的离散方法 目前最为广泛应用的为河道一汉点一河道模型，求解联立方程组的方法也有了 多样化发展。主要求解方法有：有限差分法、特征线法和有限体积法。 1 3 ，2 1 特征线法 顾元棱等人将圣维南方程组转化为相应的特征方程组，再用不等距偏心插值特 征线格式对其进行差分离散，得到内点和外点的差分方程，再在汉点处利用s t o c k e r s 条件，进而得到方程组的数值解【捌。 1 3 2 2 有限差分法 由美国国家环境保护局环境研究实验室开发的w a s p 5 ，其水动力模型程序 d y n h y d 5 采用显示格式i 硐。但显示格式是条件稳定的，其稳定性受柯朗条件的限 制，对计算时间步长要求非常严格，而且显式格式会带来新的问题，即相邻汉点之 间以及汉点与其相邻的内断面之问，流量和水位会存在一定不协调，需要采取相应 的处理技巧【2 9 1 。隐式差分格式具有稳定性好、计算精度较高、可以取较大时间步长 等优点，因而在河网计算中得到广泛的应用。对于大型环状河网，用隐式格式求解 将形成一极不规则的大型稀疏矩阵，计算量大且计算精度难以保证，给数值求解带 来极大困难。因此，河网非恒定流计算主要是针对如何减少内存占用量，如何缩短 计算时间和提高计算精度而进行的。 7 矩阵标识法在感潮河网非恒定流计算中的应用 ( 1 ) 直接解法 直接解法的基本思想是以河道断面的水力要素为基本未知量，直接求解由内断 面方程和边界方程组成的方程组。李岳生从河网矩阵的特点出发，提出了“河网隐 式格式的稀疏矩阵解法”f 3 0 j ，能够有效地节省内存，提高计算速度。 ( 2 ) 分级解法 分级解法首先由荷兰学者d r o n k e r s 提出，其基本思想是先将未知水力要素集 中到汉点上，待汉点未知量求出后，再回代求解各河段未知量。按方程组的连接形 式不同，分级解法又可分为二级解法【3 i j 、三级解法【3 2 】、四级解法【矧、汉点分组解法 m 3 5 1 等。二级解法求解所有边界条件和河段方程构成的二级连接方程组，可得到河 道首末断面水力要素，然后回代各微段方程，即可求出所有河段内部断面的未知量。 对二级解法的河段方程自相消元，可以得到以节点水位为基本未知量的三级连接方 程组，这就是目前常用的三级解法。 在三级解法的基础上，进一步从三级连接方程组分离出外边界方程和汉点能量 方程，可得到以河网内部汉点水力要素为基本未知量的四级连接方程组。 、( 3 ) 汊点分组解法 汉点分组解法是在分级解法的基础上发展起来的，其特点是根据实际工程需 要，灵活方便地将河网中的汊点分为任意多组。在计算中除了考虑连接本组汊点河 段水力要素的变化对该组汉点的水位影响外，该组汉点水位还受与本组相连的河段 的水力要素变化的影响”。 侯玉在此基础上提出了汉点分组新解法，其基本思想是利用矩阵的分块计算技 术将一般分级解法形成的原汉点水位关系应用于汉点分组，从而简化了递推过程”。 ( 4 ) 矩阵标识法 矩阵标识法【3 6 i 是在三级分组解法的基础上，根据节点水位方程系数矩阵的高稀 疏性，对非零元素进行标识。按照代码指示，将非零元素用一维数组存储，求解时 只对非零元素运算。此法在求解矩阵时完全避免了零元素的运算，可不考虑节点的 优化编码，在实际运用中，可以节省大量低效率且易出错的重复劳动。 ( 5 ) 非线性方法 徐小明将树状河网的松弛迭代法f 3 7 1 推广到环状河网，将任意复杂的河网视为一 系列的单一河道，而每一单一河道所形成的线性方程组系数都为五对角矩阵。将汇 流点各支流的流量作为旁侧入流处理，先给出其预估值，然后再用松弛迭代法来逐 矩阵标识法n - 惑潮河网非恒定流计算中的麻用 步逼近其精确值1 3 8 1 。将汉点处各支流的流量作为旁侧入流的处理方法，其在物理意 义上不甚明确。 1 3 2 3 有限体积法 文献 3 9 4 1 采用有限体积法和乘方律插值分布，取代了计算水力学中常采用的 “线性差分插值”。从试探水位出发求解动量方程得到一系列流量值，再利用水位 校正方程得到校正水位，作为新的试探水位进行计算，直至得到满足精度要求的流 场。该法有效地解决了数值计算的稳定性和收敛性问题。求解方法上采用迭代法， 较好地解决了离散方程组的非线性和相互关联性问题。但该法在汉点的处理上未充 分考虑汉点与多条河道的连接关系及汉点对河网流场的影响。且使用交错网格，在 边界及汉点处造成处理上的不便。 1 4 论文总体结构设计 由于矩阵标识法避免了矩阵中所有零元素的计算，无需进行河网的优化编码， 当概化河网有变动时无需重新编码及重新整理数据，具有编程简单、计算效率高， 可移植及通用性强的特点。而珠江三角洲河网河道纵横，中小河涌繁多，给河网的 概化带来一定的困难，在河网模拟的前期工作中常可能对概化河网进行微小的修 改。 综上所述，可知采用基于三级联解的矩阵标识法物理意义明确，且可大大提高 河网求解的效率，因此，本文试图研究矩阵标识法在珠江三角洲感潮河网计算中的 应用，采用一维圣维南方程组进行模拟，通过三级联解方法推导递推方程组，用矩 阵标识法求解只包含汉点水位的连接矩阵。避免了重复进行河网优化编码的繁琐工 作，也无需对数据文件进行重新整理和排序，大大减少了工作量和重复劳动，提高 了工作效率。 本论文共分为四章。第一章绪论，提出了本文的选题背景及研究意义，介绍 了当前国内外有关研究的现状，并概括了本文的主要内容。第二章介绍河网的计 算原理及方法，包括运用三级联解的矩阵标识法求解连接矩阵。第三章将应用上 述方法建立的模型应用于珠江三角洲感潮河网得到计算结果并进行模型验证。第 四章结语，总结、讨论本文的研究内容。 9 型型型型型型塑趔龇 第2 章河网水动力计算中的矩阵标识法解法 本章根据前述的原理与方法，重点论述河网水动力模型的建立与求解涉及控 制方程组的简化、方程组的离散和求解、初始边界条件的确定、模型的率定和验 证等问题。 2 1 河网的特性分析 河网非恒定流的数值计算，是水利和航运部门根据工程的要求经常需要避行的 工作。在现有的各种一维非恒定流数值计算方法中，隐式差分由于无条件稳定 对时间步长的选取和河段划分等方面限制较小而常应用。但在河网情况下。蓉蛐 道交叉衔接，所形成的系数矩阵是一个不规则，不对称的稀疏矩阵，需要存量较 大，且矩阵中无效零单元所占的比例也大，降低了运算速度，也相应降低丁隔许 差分方法的使用价值。 图2 一l 象征性河网 考察如图2 1 所示的象征性河网，其中：、各点为汉点。 、各点为边点。相邻两汉点之间的单一河道为河段，河段内两 1 0 矩阵标识法在感潮河网非恒定流计算中的应用 个计算断面之间的局部河段为微段。一个河段可以含有一至多个微段。微段是计 算中的基本单元，其数量决定于计算要求。如果一河网的断面总数为n 。，河段数 n ，汉点数n j ，边点数n b ；则该河网含有n s - - n ，个微段，需要求解的未知数是各 断面上的水位和流量共2 n 。个，相应需要2 n s 个独立的计算方程。 2 2 三级联解法 根据河网方程的特点，采用三级联合解法进行河网的水动力计算。即将参加 计算的方程分成微段方程和汉点连接方程，通过微段方程先求得只含有汉道断面 未知水位和流量的河段方程，并代入汉点方程以消除其中的流量或水位变量，从 而得到一组仅含有汉点流量或汉点水位的联立方程组。只要求解这一方程组，就 可以通过回代得到所有断面上的未知量。通过参考文献 3 a ，4 2 ，可推导出河网计 算方程组如下。 2 2 1 控制方程及差分格式 1 、控制方程 控制方程为圣维南方程组包括连续方程和动量方程。一维非恒定流基本方 程组为： bo _ z + 塑口 ( 2 1 ) o to x i o a + 昙( 争+ 四豢+ g 缘= 。 ( 2 - z ) 式中，z 为断面水位，q 为流量，g 为重力加速度，以为不同水位下的过水 断面面积，口为旁侧入流流量，r 为水力半径，t 7 为谢才系数，x 、f 为空间和时 问坐标。 2 、差分格式 采用p r e i s s m a n n 四点偏心隐式差分格式离散方程组： 矩阵标识法存感剃州删非恒定流计算中的应用 r l m ；u 知+ f ? ) 2 剖，口( 壁k 笪x ) + ( 1 卅( 譬)缸l 、“ 缸 雄 嚣+ f ；n f ；h 一 ； o t i m 2 a t 该格式为隐式，具有恒稳定性。 2 2 2 各级计算方程及求解步骤 1 、微段方程 将离散格式( 2 - 3 ) 式代入连续性方程( 2 1 ) 式和动量方程( 2 - 2 ) 式得： - q ；+ 1 + q m n + l + c ，z ；“+ cj z 嚣= d j e i q j “七g i q ：3 一f i z ；。+ fj z ：叠= 甲| 其中铲缶a 1 a x j ，+ 。，。q ：；1 a x j 一1 - 。o ( o j - ，一q ；) + 。，( z a 。+ z ；) b b ；+ 刚7 2 铲急一+ r 盟2 0 c 2 r 恤j ，e ，。盖一 ) ：“卫p 缸， g i = 急+ 一( 费：+ l 血，。盂+ 似) ：+ l + 互翔；+ - 血， f j = ( ) “ 矿盖恻) - 字 ( 一刎 - 了1 - 0 ( w 叫) 由曼宁公式c ：6 删基；裂 ，l伽k 删 ( 2 - 3 ) ( 2 4 ) ( 2 5 ) 矩阵标识法在感潮河嘲非恒定流计算中的应用 为书写方便，忽略上标n + l ，可把式( 2 4 ) 、( 2 5 ) 的任一微段差分方程写为： 一g + q + l + c j zj + c j z = d j( 2 - 6 ) e ，q j + g q 川一z ，+ f j ，+ l ；妒， ( 2 - 7 ) 其中。、d 、e ，、f i 、g ，、妒，均由初值计算，所以方程组为常系数线性方 程组。对一条有m 一1 个微段的河道，有2 ( m 1 + 1 ) 个未知量，可以列出2 ( m 1 ) 个方 程，加上河道两端的已知边界条件，形成封闭的代数方程组。 断面编号1 23 ii + 1 i n 一1m 单元河段编与了卜f _ _ _ _ + 了+ _ = = 吲 q 1 = ，l ( z 。) 一q l + q 2 + c 1 z i + c 1 2 2 = d l 臣q l + g 1 q 2 一f l z l + e z 2 = 伊1 一q 2 + q ，+ c 2 2 2 + c 2 2 3 皇d 2 e 2 q 2 + g 2 q 3 一f 2 2 2 + e z 3 = 妒2 一q ，+ q ，+ 1 + c j z ，+ c i z ，+ 1 ；d ，( 2 - 8 ) e i q i g i q i n f | z i + f i z | n = 甲i q 。一1 + q ，+ c m 一1 z 。一l + c 。一1 z 。= o m e 。一1 q l + g 。一1 q 。一f 二一1 z ，一1 + f i i z 。昌妒。一l 瓯= 厂( z 。) 这一类方程，每一方程只含有四个非零系数，并且排列整齐，未知量的编号 前后相连，在总的系数矩阵中形成一个紧邻于主对角线的阶状分布。在个河网 中，该方程的总数是2 ( n , - n 0 个，是计算方程中最多的一类。 2 、河段方程 对河网中的任一单独河道，其微段方程都是前后相连有序的，可对其自相消 元得到只含有河段首、末断面变量的河段方程组，即只含有节点变量： q 1 = a l + 芦1 2 1 + 6 1 z 。( 2 - 9 ) q ，= 日。+ 叩。z 。+ r 。z l ( 2 一l o ) 其中z 。为首节点水位，z 。为末节点水位。 3 、节点连接方程 在节点上，水流总满足以下两个方面的条件： 矩阵标识法在感潮河嘲非恒定流计算中的应用 ( 1 ) 流量衔接条件：进出每一节点的流量必须与该节点内实际水量的增减率 相平衡 q j 。票( 2 - 1 1 ) 式中，f 表示节点中各节点断面的编号，q f 表示通过i 断面进入节点的流量。 q 为节点的蓄水量。如将该点概化成一个几何点，则q ：o 。 ( 2 ) 动力衔接条件：如果各断面的过水面积相差悬殊，流速有较明显的差别， 当略去节点的局部损耗时，由伯努力( b e m o u l i ) 方程得： z 1 + 芒。z 2 + 譬 ( 2 - 1 2 ) 2 9。冶 、 如果节点可以概化成一个几何点，出入各节点的水流平缓，不存在水位突变 的情况，则各节点断面的水位应相等，且等于该节点的平均水位，即 z 。一z ：s 一z( 2 1 3 ) 在更一般的条件下，如汉点处设有闸、堰等建筑物，汉点的动力衔接条件总 可以按具体条件概化成一般形式： a h i + 6 q f + c h j + d q 。= e ( 2 1 4 ) 在每一个汊点上，如果含有n d 个汉道，则可以列出两两汉道之间独立的动力 衔接方程n d 一1 个，流量衔接方程1 个，合计n d 个，即等于汉道的总数。而每一河 道总有两个端点，不属于边点便属于汉点，因此一个河网内汉点连接方程的总数 是2 n ，。n b 个。这类方程中的非零系数虽然也不很多，最多不超过2 n d 个，但其在 系数矩阵中的位置取决于汉道断面的编号；而一汉点上的断面编号一般说来不仅 不能相连成序，并且很不规则，造成矩阵中系数的排列也很散乱。因此在一般求 解方法上，汉点编码存在优化编码的问题。 4 、边点方程 边界条件包括水位边界、流量边界以及水位流量关系3 种，对于一个大型河 网的计算，比较理想的边界条件是既有流量边界，又有水位边界条件。尤其是对 于受潮汐影响的河网，仅有水位边界条件是不足以得出准确的模拟结果的。 边点方程式有以下几种情况： 1 4 矩阵标识法在感潮河蹦非恒定流计算中的应用 a 。= q 。( f ) z ，= z 1 8 ) a ，= f ( z 。) 相应的控制方程总可以概化成如下形式： 口z f + 6 q f = c ( 2 - 1 5 ) 其中。、b 、c 为按实际条件导得的系数和右端项。这一类方程的非零系数最多只 有两个，分别位于主对角线及紧邻该线( 前或后) 位置上，排列也是比较整齐的。 一个河网中的边点方程总数n b 个。 微段方程、节点连接方程和边点方程构成了求解的闭合方程组，方程总数为 2 m 个，等于未知数总量，因此可以联合求解。 5 、求解步骤 ( 1 ) 就各微段建立圣维南方程组并经过隐式差分近似得微段方程组： ( 2 ) 将每一河段内的微段方程依次消元得到河段方程； ( 3 ) 将河段方程改写成水位或流量的隐函数，代入边点方程和汉点连接方程 得到以流量或水位为变量的连接矩阵； ( 4 ) 求解连接矩阵得到各汉点上各断面的流量或水位； ( 5 ) 回代河段方程得汊点各断面上的水位或流量： ( 6 ) 回代微段方程得所有各个断面上的水位和流量。 整个求解过程，可以用框图表示为：、 1 输入原始数据2 形成连接矩阵 指示信息数组 3 逐时步建立微 段方程组 4 建立河段方程l 一5 确定边界方程l 一+ f6 建立连接矩阵 7 连接矩阵消元 得水位 8 阿代河段方程 得汉点流量 图2 1 三级联解计算框图 1 5 9 同代微段方程 得全部未知量 1 0 结束 矩阵标识法在感潮河耐非恒定流汁算中的应用 框图以河段方程中的流量为隐函数，形成的连接矩阵以水位为未知量。也可 以水位为隐函数，建立流量连接矩阵，其求解步骤完全相同。 2 2 3 求解过程中的关键问题 2 2 3 1 河段方程的推导 对于河网，内河道和外河道的计算方法是不同的。对外河道，河道的一端是 边界已知的外节点，而另一端是水力要素未知的内节点。内节点的变量必须通过 河网的联解才能得出。为了减小方程组的阶数，通常以内节点构成基本河网。 1 、外河道的追赶方程 外河道的求解是，通过边界条件确定首断面的递推关系，用追赶的方法求解 各断面的水位与流量递推关系，得到汇入基本河网的末断面与内节点水位关系。 当内节点水位求出后，即可回代求出各断面的水位和流量。当边界条件不同时追 赶方程也不一样。 当上游边界为水位条件已知时，有追赶方程： z ，一只一k q l( 只一z ，( f ) ，k = o )( 2 - 1 6 ) 代入( 2 _ 8 ) 的第一单元微段得： 方程中有3 个未知数q l ，q 2 和z ：，为了解出幺、z ：，需要用q ：表示q 1 、z ：。 设有如下追赶方程： 腭- 妥一罂( 2 - 1 8 ) l z ：= b 一心q ： 式中系数s ：，疋，只，心可由( 2 1 7 ) 式的系数唯一确定，将( 2 1 8 ) 式代入( 2 8 ) 式的第二单元微段，得： j 罟：2 1 s 一罂 ( 2 - 1 9 ) l z ，= 只一屹q 3 、 以此类推，对第i - 1 个单元微段有： 1 6 乃 q q = = 乞乙g 卜卜g麓假e q一： + q 她崛蜴p 矩阵标识法在感潮河嗣非恒定流计算中的应用 l q 一1 = 6 ，一，世j i z j = p 一匕q ， 将( 2 2 0 ) 式代入( 2 - 8 ) 的第j 单元微段方程，得： f - q j + q 州+ c ，( p j 一9 ，) + c ，z 。d j 1 e j q ，+ g q 一f j ( 弓一v , o j ) + f j z j “2 妒， 以q 。为自由变量，求解得： 嚣篇臻+ 1 ( 二，m - d 舯s ? 等等一。= 丽c i g _ ，- f j ， 一竽一。；半； y l = d i - c i p l ，y 2 一甲i + f i p ! ，y j ；i + c i v f ，y 4 。e i 七f y i 。 当上游边界为流量条件已知时，同理可以推得追赶方程如下： q 1 一只一k z ，( 只= q 1 ( f ) ，k = o ) f z ，= s t + 1 0 + t z 川 i q ；鼻r _ 十l z 鼽。= 丽a y 3 - 1 4 一，= 鬻， 只+ 。= y 3 一k s ，“= c ，一k t “； k 一+ c j ，砭；f j + e ，y 3 一d j + 0 ，y 42 9 ，一e j e j ( 2 2 0 ) ( 2 2 1 ) ( 2 - 2 2 ) ( 2 2 3 ) ( 2 2 4 ) 2 、内河道的追赶方程 对于内河道，首、末断面的边界条件均为未知，由于没有端点边界条件可供 利用，因此不能适用外河道的追赶方程。设内河道中任一河道首断面号为1 ，末断 面号为m 。 用逆推法( 从河道的最末一单元微段方程开始往前推导) 可以把内河道中任 一断面的流量表达成本断面水位和末断面水位的线性函数。由( 2 8 ) 式最后单元微 段方程消去q 。，解出妇 棚一，得： 1 7 矩阵标识法在感潮河网非恒定沉计算中的应用 q _ 一l = 口，一1 + 卢。一l z 。一1 + 。一1 z 。，( 2 - 2 5 ) 其中，a 。一，= ! ! 警，卢。一。= ! j 筹， 。一t = 兰：糌。 设己推出如下式子： q + 1 一口j + 1 + “z ，+ l + 言+ i z 。， ( 2 - 2 6 ) 将之代入( 2 8 ) 式第f 单元微段方程有： 一q j + ( a + l + 卢，+ l z ，+ 1 + 亭+ 1 z m ) + c ，z ， l c j z j + i = d i ( 2 - 2 7 ) i e ，a j + g ( a j + l + 芦j “z “+ 毒，+ l z m ) 一e , z + zj + l = 妒， 解出o ；得： a i = a l + p i z i + l i z 。 鼽叩迎气宥竽型， ( 2 2 8 ) 铲警， ，考川吼一x g ，) 钉4 1 霭万 k = c + 卢j “，k 置g ，卢，+ 1 + f j ， ( j = m 一2 ，m 一3 ，d 同理，用顺推法( 即从第一一单元微段方程开始往后推导) 可把内河道中任 断面的流量表达成本断面水位与河道首断面水位的