二维水沙数学模型的建立和应用ppt课件.ppt

天津大学硕士学位论文平面二维水沙数学模型的研究与应用 1 摘要二十世纪70年代以来 随着电子计算机的普及与计算机性能的不断提高 人们在研究河流泥沙问题时 越来越多的采用了数学模型 它以数值方法和计算机技术为手段 通过对河流的数值模拟计算 解决河流工程所关心的问题 对二维水流泥沙数学模型的探索和研究已在水利 环境 港口航道 防洪等工程领域中取得了一定成果 但水流泥沙运动机理仍是工程问题中的重要研究课题 本文结合河北省水利水电勘察设计研究院 河北省庙宫水库汛期调度运用方式 项目 研究了平面二维水流泥沙数学模型的理论 方法和实际应用问题 主要工作如下 1 系统地推导了二维非恒定流 二维泥沙连续方程 应用有限体积法思想 采用分类简化离散方程的方法 降低了计算难度 节省了计算时间 为建立平面二维水沙数学模型奠定了基础 2 在已有的泥沙研究成果的基础上 提出并推导了新的泥沙连续方程形式 结合适合庙宫水库情况的爱因斯坦水流挟沙公式 建立了平面二维水沙数学模型 2 3 根据离散方程 编制Fortran计算程序 对庙宫水库库区内的水流和泥沙冲淤进行计算并与已有的实测资料对比 吻合效果较好 说明模型可以用于该水库的河道变形模拟 对上游来水 来沙情况和水库排沙控制水位序列13年水沙组合过程进行计算模拟 给出了全河段的流场分布 含沙量分布 河床冲淤演变分布的计算模拟成果 4 通过综合对比 定量分析了不同计算组合下床面冲淤演变情况 计算出排沙比 提出了水库现状泄流条件下最优清淤方案 研究成果不仅为理论研究提供了新思路 而且也为实际工程提供了可行的数值模拟方法 具有一定的理论创新和应用价值 关键词 平面二维水流泥沙数学模型 庙宫水库 泥沙连续方程 冲淤演变 3 目录 第一章绪论1 1泥沙数学模型发展简况1 2水沙数学模型基本方程及计算方法1 3数值计算方法1 4水沙基本理论存在的问题1 5淤积上延现象第二章河道复杂边界的处理技术2 1斜对角笛卡尔方法的基本原理2 2斜对角笛卡尔方法的优点2 3边界条件及边界网格单元类型处理2 4本章小节 第三章平面二维水流泥沙数学模型的建立3 1二维水沙数学模型的基本方程3 2初值条件及边界条件的处理3 3网格划分3 4水流挟沙能力第四章平面二维水沙数学模型的数值模拟4 1水流控制方程的离散4 2泥沙连续方程的离散4 3河床变形方程的离散4 4模型求解过程 4 第五章水流泥沙数学模型验证及应用5 1模型的基本特征5 2网格划分5 3模型验证5 4模型预报5 5排沙量及冲淤分析第六章工作与展望6 1工作6 2展望 5 目前研究河流泥沙及河床演变问题 主要有三种手段 即原型资料分析法 物理模型方法 数学模型模拟方法 这三种方法各有千秋 视具体工程问题 可独立运用或联合运用 物理模型方法是将原型的尺寸以及各种水力要素在遵循相似原则的基础上按一定比尺缩小为模型 物理模型方法具有直观 准确 易于操作等特点 是早期研究河床演变的主要方法 近年来随着测量技术的进步 物理模型的精度得到了很大的提高 但是与数学模型相比 物理模型耗资大 对外界条件要求高 且完成实验方案的周期长 河流泥沙数学模型能够模拟水流泥沙运动过程 河床冲淤量 河床形态及河床泥沙组成的变化 它的研究和应用可以追溯到上世纪50年代 由于受到当时计算工具的限制 这种计算仅限于一些简单的状况 而且精度较差 工作量大 因此 未能得到普遍推广 二十世纪70年代以来 随着电子计算机的普及与计算机性能的不断提高 这些困难已逐步克服 数学模型的优越性越来越突出 经过近二 三十年的发展 水流泥沙数学模型经历了由一维到二维 到二 第一章绪论1 1泥沙数学模型发展概括 6 维嵌套 到三维 到三维嵌套 由原先在概化水文 泥沙及河床条件下数学模型的建立 率定和验证 发展到与物理模型相互配合共同回答工程水流和泥沙问题 复合模型 的阶段 1 2水沙数学模型1 水流控制方程2 泥沙控制方程 7 8 1 3数值计算方法概述1 3 1数值离散方法 描述河流水沙运动的控制方程一般是非线性的 对于这些控制方程 一般很难求得解析解 只能通过离散控制方程 求解代数方程组来得到近似数值离散解 数值计算方法的精度及速度依赖于控制方程的离散方法 代数方程组的求解方法 网格的划分及边界条件的处理等 目前常用的数值计算方法有 特征线法 有限差分法 有限元法 有限体积法 有限分析法等 1 特征线法特征线法是河流模拟数值计算在计算机尚未普遍之时所采用的一种方法 主要采用手工计算 其最初思路是在平面上绘制特征线 其优点是能反映问题中信息沿特征传播的性质 算法符合水流运动的物理机制 稳定性好 计算精度高 但对于求解周期短 变化急剧的问题计算效果较差 因此目前很少直接用于数值计算 但是其原理仍很重要 经常用于作为了解其他数值方法的基础 2 有限差分方法有限差分方法 FDM 是计算机数值模拟最早采用的方法 至今仍被广泛运用 FDM以Taylor级数展开等方法 把控制方程中的微商用差商代替进行离散 从而建立代数方程组来求解该方法数学概念直观 表达简单 其解的存在性 收敛性和稳定性早已有较完善的研究成果 是比较成熟的数值方法 目前应用广泛 9 3 有限元方法 有限元法是根据极值原理 变分或加权余量法 将问题的控制微分方程化为控制所有单元的有限元方程 把总体的极值作为各单元极值之和 即将局部单元总体合成 形成嵌入了指定边界条件的代数方程组 求解该方程组就得到各节点上待求的函数值 有限元的求解方法常见的有直接法 变分法 加权余量法和能量平衡法等 其优点是网格划分灵活 拟合复杂河岸边界容易 网格节点可局部加密 稳定性好 精度高 适合于几何 物理条件复杂的问题 对于隐式FEM 其精度较高 但数学推导繁杂 计算量和储存量较大 而且在误差估计 收敛性和稳定性等方面的理论研究与有限差分法相比还显得不够成熟和完善 尤其在多维计算中 由于有限元法贮存量大 会直接影响计算速度 10 4 有限体积法 有限体积法是计算出通过每个控制体边界沿法向输入 出 的流量和动量通量后 对每个控制体分别进行水量和动量平衡计算 便得到计算时段末各控制体平均水深和流速 因此 有限体积法正是对于推导原始微分方程所用控制体的回归 因为控制体间界面输运的通量 对相邻控制体来说大小相等 方向相反 故对整个计算域而言 沿所有内部边界的通量相互抵消 对由一个或多个控制体组成的任意区域 以至整个计算区域 都严格满足物理守恒律 不存在守恒误差 并且能正确计算间断 由于采用守恒型的微分方程并对每一计算单元进行质量和动量守恒形式的离散 使得微分方程包含的守恒性质在每一个控制容积上都得到满足 若保持各单元界面两侧相邻控制体的计算输运通量相等 那么整个计算区域上都能保持守恒 11 在数值模拟的过程中 虽然采用的离散求解方法不同 但都有相同的特点 即首先把计算区域划分成许多控制体或网格 然后在这些小控制体上把微分方程离散成代数方程 再把小控制体上的代数方程汇合成总体代数方程组 最后在一定的初始条件下求解此方程组 从而求得计算区域内各节点的物理量 所以数值模拟的正确性和精度取决于网格的划分 方程的离散 初始值条件 代数方程组的求解以及所建模型的物理理论依据是否正确合理等几个因素 各种方法均有其自身的优点和适应性 在实际计算时选择什么数值方法应根据所研究问题的特点和计算精度要求 以及研究者的习惯而定 12 1 3 2边界处理技术 1 不规则边界的处理在水利工程中 水流的边界往往是很复杂的 而在数学模型计算中如何处理不规则的复杂边界将直接影响到整个模型计算结果的可靠性 精度和计算时间 是数值模拟的一大难点 选择复杂边界的处理方法与算法的选择有着同等重要性 常用的不规则复杂边界处理方法主要有 美国学者Thompson 16 提出的贴体坐标法 非结构性网格以及多重网格法 贴体坐标法因其能使计算域边界上的坐标线与边界线密切贴合 而获得了广泛的研究与应用 其主要优点是用曲线坐标逼近复杂边界 在边界模拟上具有较高的精度 但该方法在较复杂的边界上不易形成正交贴体网格 如果是非正交曲线坐标系 则会产生大量的交叉导数项 给求解带来很大的困难 非结构性网格使用三角形和四边形计算单元 使得在复杂区域形成计算网格更加灵活和方便 但非结构性网格往往要建立一套标明网格节点之间关系的复杂系统 因此在编程计算处理时比较麻烦 特别是针对自由面流动和干湿边界变化时 这种困难更加突出 多重网格法则由于网格系统之间的相对独立性 难以做到对系统守恒性的保证 另外对于复杂边界的处理问题 俞国良于1987年提出了镜像法 周建军于1988年提出了边界对称点法 这两种方法都是根据对称原理 假设流场外有一定的虚拟流动 从而导出计算域外节点各物理量的计算式 该方法简单易行 提高了不规则边界上运用差分方法的计算精度 但由于域内外对称点一般不在节点上 仍需一定插值过程 会带来一定误差 13 2 动边界的处理 由于河床冲淤变化和水位的波动 在河岸及洲滩等区域都会遇到动边界问题 动边界数值模拟的困难主要在于沿动边界法向流动不同于明渠均匀流 其非恒定和非均匀性强 常用的曼宁摩阻公式等在形式上难以套用 再者是因为水深很小 对离散格式的要求很高 要求数值解不产生假振 保证水深总是大于零 常见的处理动边界方法有两种 一种是对实际边界的数值模拟 即时刻追踪动边界的准确位置 另一种方法是 把整个计算区域都按照常进行计算 在计算过程中取一定的方法对无水区域进行处理 但是由于第一种方法在编制程序中处理较复杂 因此目前经常采用的方法是第二种方法 在计算过程中不追踪动边界的准确位置 把整个计算域都参与运算 对于干河床部分采用一些处理技巧进行处理 比如何少苓提出的 窄缝法 其思想是在床面上人为设置窄缝 由于缝足够深 总是有水流动 而缝很窄 当河床露出水面时由窄缝引起的水量守恒误差很小 此外 还有周建军提出的渗透边界法 程文辉 王船海等运用的 冻结 技巧 14 1 4水沙基本理论中存在的问题 1 阻力问题阻力计算是阻力系数 糙率及涡流粘性系数 的确定问题 处理阻力的方法有两种 一种是采用水力半径分割法或能坡分割法划分阻力单元 分别计算各单元的阻力 然后计算总阻力 另一种方法是根据实测资料直接计算总阻力 以上两种方法都只能用于一维水流糙率的计算 多数的二维数学模型均是直接采用一维阻力系数进行计算 这种处理方式显然比较粗略 李义天在考察一维情况下断面综合糙率沿河宽变化分布的情况下 提出了二维糙率计算公式 实际应用表明 该方法比直接采用一维方法确定的糙率更合理 能比较好地反映实际情况 但由于因素比较复杂 由此得到的推论在定量上尚未能被实测资料所证实 2 水流挟沙力问题水流挟沙力的计算一般有两种方法 一种是将含沙量沿垂线分布公式积分求垂线平均含沙量 但河底泥沙含沙量很难确定 不同的研究学者给出了不同的计算方法 所得的计算结果与实际情况差别较大 另一种是半经验半理论方法 这样的公式很多 公式中包含两类参数 一类是反映水力条件的参数 如流速水深等 另一类是反映泥沙特性的参数 如沉速 泥沙粒径等 这种方法是建立在实测资料的基础上得 因而结果比较接近实际情况 得到了广泛的应用 现在挟沙力的研究成果多基于一维问题 对二维问题研究的较少 而大多数的二维问题也是直接采用一维挟沙力公式进行计算的 但是李义天通过经验比较 发现二维挟沙力与一维挟沙力有很大的差异 15 3 床沙悬移质沿程变化的问题 目前在数学模型计算中存在的尚未很好地解决的问题是床沙质及床沙级配的沿程变化问题 对于冲淤甚小的处于自然状态下的河流 在许多情况下床沙及其级配在冲淤过程中可以假定基本不变 而冲泄质级配则视来水来沙条件而定 对于冲淤幅度甚大的河流 对床沙质及床沙级配的沿程变化不能不进行考虑 现在大多数模型将推移质和悬移质分开 将悬移质中床沙质与冲泄质分开的做法基本相同 至于将床沙质进一步分级来求不同粒径的冲淤变化 一般都比较粗糙 韩其为 23 提出了一套考虑泥沙级配变化的计算模式 其基本思想是不区分床沙质与冲泄质 假定水流挟沙力的级配和实际输移的级配一致 而悬移质和床沙质级配在每一时段内的变化都看成本时段内河段冲淤的直接后果 计算断面每个粒径组的含沙量须反复试算才能确定 该方法对泥沙输移的物理考虑比较细致但计算过于复杂 16 1 5淤积上延现象 当水流进入到水库回水区后 由于水流条件的改变 流速减小 泥沙落淤 而落淤的结果反过来又促成水流条件的进一步改变 使得回水曲线在淤积区及其上下游一定范围内普遍抬高 为泥沙继续在原淤积区落淤并同时向上下游发展提供条件 这种淤积自回水末端 通常指正常蓄水位回水末端 向上游发展的现象 称为水库淤积的上延现象 也就是通常所说的 翘尾巴 现象 按淤积中剖面形态的类型 可以分为三角洲淤积 带状淤积 和锥体淤积见 17 图1 1官厅水库淤积纵剖面 图1 2丰满水库淤积纵剖面 18 图1 3陕西黑松林水库淤积纵剖面 19 前的3 74 变缓为1 8 淤积末端距大坝13 7km 上延系数1 55 随着淤积末端不断上延 上游河床逐年抬高 距大坝9 6km的四合永公路桥处河床淤高3m多 使河床由原来的地面以下到现在与两岸持平 甚至高出两岸耕地及村基 产生严重浸没 浸没区内数千亩耕地盐渍化 沼泽化 不能正常耕种 近千户居民住房地基沉陷 墙垣断裂甚至倒塌 给群众生产带来严重损失 以至出现影响社会稳定的问题 并严重威胁四合永镇的安全 受到省 部领导和有关部门的高度重视 为延缓庙宫水库淤积浸没的发展 自1987年至1992年连续6年进行了空库排沙原型实验 试验中对不同水沙条件 空库排沙时间 不同排沙方式等进行了全面的试验 总结分析了水库现有泄流排沙设施条件下的减缓水库淤积浸没的可能性和运用效果 1993年由水利部天津院完成的庙宫水库空库排沙试验物理模型研究和同年由天津大学泥沙研究室完成的庙宫水库排沙洞前期工程泥沙数学模型初步研究 为当时增设新排沙洞位置 洞径和底高程的确定提供了理论依据 近年来随着水环境关注意识的提高和水资源高效利用思路的确立 如何在水库保水和水库排沙间选择优化的水库运行模式是需要进行研究和解决的重要问题 根据庙宫水库设计正常蓄水位778m 大沽高程 以下不加说明处均采用大沽高程 和库区死水位768m的要求 优选即能泄水排沙又能保水兴利的调水调沙方案是本次研究的主要目的 研究过程拟采用库区限蓄水位768m 770m 772m 774m 776m 778m六级水位标准 来水来沙过程采用1958 1959 1971 1975 1985 1989和1990的四个系列三种组合来模拟不同年份水沙过程 在18个计算系列中进行方案优选 20 第二章河道复杂边界的处理技术 自然界中 水流的边界往往是十分复杂的 而在水流数学模型计算中如何处理不规则的复杂边界将直接影响到整个模型计算结果的可靠性 精度和计算时间 是数值模拟的一大难点 复杂边界的处理方法与算法的选择有着同等的重要性 因此在研究江河湖泊的水流数值计算过程中 必须要对边界作出一个适当的处理方法 2 1斜对角笛卡尔方法的基本原理斜对角笛卡尔方法基本原则就是按就近逼近原则 把距实际水边线较近的网格节点连线以拟合实际水边线 生成网格时 根据边界是否固定 可采取不同的判定方法 当边界为固定边界时 可把固定边界离散为若干离散节点 确定各离散节点坐标 按距离最小原则选定固定边界上各离散节点的逼近网格节点 然后依序连接各网格节点 即形成拟合边界 对于动边界 由于边界线变化不定 则可借助网格节点间水深的干湿变化和水深值 参照相应网格节点间的地形高差 即可判定距离实际水边线较近的网格节点 合并纵横向的拟合网格节点 即可形成拟合瞬时水边线 21 图2 1斜对角笛卡尔法边界拟合示例 22 2 2斜对角笛卡尔方法的优点 1 边界拟合效果良好 更接近实际边界 由图2 2可以看出 计算水边线与实际水边线拟合良好 其边界长度误差一般小于6 仅为锯齿法长度误差的 且随网格间距的缩小而更为精细 图2 3为按锯齿法 也以距离最近为原则 确定的计算水边线 对比图2 2 图2 3可知斜对角笛卡尔法所拟合的水边线比锯齿法拟合的水边线更为贴近实际水边线 事实上 斜对角法就是在锯齿法基础上通过斜对角连线添加一个三角形网格单元 使边界形态更接近真实形态 且避免了锯齿状边界的产生 使流场计算更精确 2 处理简单 附加计算量小在矩形网格水域边界上运用斜对角笛卡尔法 计算时 只须对三角形边界网格单元进行处理 而三角形边界网格单元若把斜边中点视为一角点 则仍可视为四边形网格 可用统一方程求解 与内部计算网格无异 并不像贴体坐标拟合一样须对基本方程进行坐标变换而增添许多交叉导数项 因此 用斜对角笛卡尔法进行边界处理 不仅保留了离散简便 计算效率高的优点 同时又使边界拟合更加精细 避免了锯齿状边界对流场计算的影响 23 图2 2斜对角笛卡尔法边界拟合图2 3锯齿法边界拟合 24 2 3边界条件及边界网格单元类型处理 根据斜对角笛卡尔法的逼近原则 在平面二维计算中拟合水边线上形成的计算边界网格单元存在八种基本类型 见下图 图2 4八种网格单元基本类型 25 图2 5网格类型及差分点分布 26 斜对角笛卡尔方法不仅可以应用于有限差分法中 也可应用于有限体积法等其它离散格式中 它对复杂边界条件和复杂水流条件的适应性和计算精确性将使复杂边界的处理变得简单高效而且准确 为复杂边界问题和近岸工程问题的计算提供了一种新思路 27 第三章平面二维水流泥沙数学模型的建立 天然河道水沙数学运动一般都属于三维运动 运动要素既沿程变化 又沿水深和河宽发生变化 28 由于三维水流运动比较复杂 难于计算 常将运动要素沿水深方向平均 把三维问题转化为平面二维问题 平面二维非恒定流方程 也就是浅水方程 它描述的对象是符合以下条件的均匀流体 1 具有自由表面 2 以重力为主要的作用力 以水流与固定边界之间及水流内部的摩阻力为主要耗散力 有时还存在水面气压场风应力以及地球柯氏力等的作用 3 水平流速沿垂线近似均匀分布 不考虑实际存在的对数或指数等形式的垂线流速分布 4 水平运动尺度远大于垂直运动尺度 垂向流速及垂向加速度可以忽略 从而水压强接近于静水压强分布 28 在实际问题中 通常可以作为非恒定流处理的具有自由表面的实际水流 有如下情况 1 水深相对较浅 可以根据水深h和波长L之比区别 当 0 4时称为浅水长波 2 水底坡度较缓 设底坡倾角为 判断缓坡的条件是 此时 可以忽略底坡引起的垂直速度和垂直环流 也不必考虑垂直加速度及由此产生的动水压力 3 水面渐变且坡度较缓 既使满足以上两个条件 在特定的情况下 由于受到天然或人工控制 有时水面比降仍很大 水深沿程变化很大的急变流 水压力不能用静水压强分布来近似 4 无明显垂直环流 由于外力改变或地质环境的改变在垂直平面内产生二次流 会给浅水模型带来误差 需要加以改正 此时方程不能沿垂向积分 要按不可压流来处理 密度沿垂向变化较大 改变了静压分布 需要作为密度流处理 因此 通常的浅水流动指在重力作用下密度均匀 具有自由表面 流动近似水平的长波传播现象 相比一维水沙数学模型 平面二维水流泥沙数学模型对沿水深平均的水沙各运动要素 垂线平均流速 悬移质含沙浓度 推移质输沙率 以及河床变形等随流程沿河宽两个方向的变化情况均可以进行比较细致的模拟 29 1 水流控制方程 3 1二维水沙数学模型的基本方程 30 2 泥沙控制方程 31 3 关于悬沙连续方程的一些探讨泥沙在水流中分布就是悬沙的含沙量 而分布的变化即为悬沙含沙量的变化 它是可以用泥沙连续方程来描述的 悬沙含沙量的变化不外乎是由对流项 上游输沙 和本地的紊动扩散项引起的 由此可以把对流项 紊动项引起的悬沙含沙量变化用两个连续方程来表示 主要思想如下 由三维泥沙扩散方程 其中 是单元体内的含沙量 可以这样认为单元体内含沙量是由上游来沙和本地床面紊动扬沙引起的 我们把上游来沙影响用来表示 称为输沙含沙量 本地床面紊动扬沙用来表示 称为扬沙含沙量 这样 可以把单元体总的含沙量用如下式子表示 32 分别根据时均 脉动的概念 可以得到两个方程 可见 输沙含沙量和扬沙含沙量都是满足三维泥沙扩散方程的 三维方程在物理意义描述上是十分准确的 但是在实际应用中 因为计算条件在河道发生冲淤变化时 悬移质底层及床面扬沙含沙量的变化是引起床面冲淤的主要因素 而输沙含沙量只是上游来沙在本地的一种扩散效应 在当时刻冲淤变化过程中并不会起到举足轻重的作用 我们把含沙量区分为输沙含沙量与扬沙含沙量也是为了在计算冲淤变化过程中体现这种关系 以使结果更符合实际情况 33 对输沙含沙量方程和扬沙含沙量二维化 34 35 36 如上所述 悬移质中造床质所占比例是要受条件限制的 同样的 床面泥沙能否起动 以及起动后是以何种形态运动 这在模型计算中 也是要予以考虑的 Shield 1936 对各种泥沙颗粒进行了临界起动实验 实测得到了无量纲临界起动剪切应力与颗粒雷诺数的关系曲线 37 38 39 40 41 42 结果表示具有波动性 表明床面有可能出现沙纹或沙波的变化形式 同时床面的变化形式与扬沙含沙量的变化相关 图中表示的床面比较规则 但对复杂的水动力条件和床面初始形态 床面变化具有拟周期性 同时变化程度非常复杂不规则 由于泥沙问题的复杂性 目前对泥沙运动机理还没有完全认识清楚 现有的各家对基本问题的认识及处理方法也不尽相同 甚至相差甚远 在20世纪70年代开始出现了大量根据实测资料 利用计算机进行回归方法计算得出的阻力和输沙力经验公式 由于有大量的实验及天然资料做基础 各种变量的覆盖范围大 这样得出的一些公式精度很高 甚至超出了许多理论公式 其原因在于现有的理论体系往往不能够把所有的泥沙的输移过程完整地 正确地反映出来 本节内容从泥沙扩散的物理角度出发 基于现有的泥沙理论推导出两个可分别计算输沙和扬沙的公式 在理论扩展 公式推导上做了一新的改进 从公式结构上看与现有悬移质泥沙连续方程结构相同 如将输沙连续方程与扬沙连续方程相叠加得到的全沙连续方程与现有悬移质泥沙连续方程具有完全相同的形式 说明方程推导正确 只是在划分输移泥沙和扬动泥沙标准上可能还存在进一步研究的空间 43 44 3 3网格划分网格的确定对模拟实际流动的可靠性同样重要 建网基本要求是 符合流动特点 易于建立 比较光滑和规则 满足精度和计算的稳定要求 便于组成节约 高效的数据结构 必要时可随时依解的梯度作适应性调整 本文为了更好的模拟计算域内的实际地形 地貌 采用的是无结构有限体积网格划分计算区域 使模型地形 地物与实际情况基本一致 一般而言 无结构网格划分采用三角形网格 因其有较为成熟的自动生成技术 故自动化程度较高 但是针对本文设计的二维河道情况 采用的是不规则四边形组成的网格 主要考虑到 一方面若完全采用三角形网格作均匀划分 反而会难以客观反映小水量时水流沿主槽流动的特点 有可能造成计算结果失真 另一方面 网格呈四边形格式 局部辅以三角形网格 这样不但可以很好地拟合边界形状 而且格子数及其边数 要比完全采用三角形网格时少 同时实际空间步长也可以取大些 45 46 3 4水流挟沙能力 1 悬移质输沙率公式的选取本文中采用的积分求输沙率的方法 选用了爱因斯坦挟沙力理论 47 48 49 50 51 第四章平面二维水流泥沙数学模型的数值模拟 本文采用有限体积法离散方程 它将计算区域划分为网格 并使每个网格点周围有一个互不重复的控制体 将待解微分方程对每一个控制体积分 从而得到一组离散方程 求解差分方程组 其求解就作为微分方程定解问题的近似解 4 1水流控制方程的离散 52 动量方程的求解 采用分类简化处理的方法 把它概括为地面型通道 河网型通道和特殊通道来分类处理 假定同一时段内同一网格水位变化不大 同时由于平面浅水中对流项作用较小 在计算中略去不计 此外 有根据堰流流量公式计算的连续堤或缺口通道 53 1 河道型通道 即通道两侧网格均为河道型网格 动量方程中保留局地加速度项 重力项和阻力项 利用差分方法离散得到河道型通道的动量离散方程 54 2 地面型通道 即通道两侧单元为陆地地面 且通道上没有堤防等阻水建筑物 考虑到库区内的地形变化 地面洪水演进主要受到重力和阻力的作用 利用差分方法离散得到地面型通道的动量离散方程3 对于高于地面的阻水建筑物 如堤防 铁路 公路等 可以概化成连续堤或缺口堤通道 其流量采用宽顶堰溢流公式来计算 离散后得到泥沙连续方程离散方法与水流连续方程离散方法相似 此处不再赘述 55 4 2泥沙控制方程的离散泥沙连续方程离散方法与水流连续方程离散方法相似 此处不再赘述 56 57 58 第五章水流泥沙数学模型验证及应用 本模型的建立主要用于庙宫水利枢纽库区泥沙变化的研究 选择了库区河道为模拟对象 在分析库区河道水流泥沙运动特点的基础上 建立了二维水流泥沙数学模型 并用实测水文资料对模型进行验证分析 在此基础之上 进行了水沙组合的计算 59 图1 1庙宫水库位置示意图 60 61 5 1模型基本特征年际间输沙量变化悬殊 庙宫水库入库最大沙量1360万吨 1959年 最小沙量61万吨 1981年 两者相差22 3倍 蚁蚂吐河下河南站1959年输沙量为1080万吨 1966年输沙量仅32万吨 变幅达33 8倍 水沙特点表现为水沙过程的不完全对应一致 有时沙峰先于洪峰 有时反之 有时大水大沙 小水小沙 也有时大水小沙 小水大沙 主要是由于暴雨强度 过程和中心位置的差异以及各支流水土流失程度不同等原因造成的 62 5 2网格划分庙宫水库数学模型自水库大坝向上游延伸 在四合永镇北河道岔口分为两支 一支沿伊逊河向西北至掌字水文站 另一支沿不澄河向东北至边墙山水文站 模型主体南北走向 全长约12 6km 模型最宽河段约1km 模型内包含27个实测断面 公路桥 铁路桥各一座 模型网格划分以实测断面为基础 采用无结构有限体积方法划分网格 考虑主河槽走向 划分后得到断面194个 单元3451个 结点3626个 通道7076个 断面平均间距为65m 断面内空间步长为50m 单元平均面积为2124m2 模型覆盖水域面积约为8km2 模型内实测断面分布见图5 3 模型划分网格见图5 4 网格单元编号 通道编号 节点号的分布形式见图5 5 63 64 65 66 5 3模型验证1 水流模型验证验证水流模型上游边界采用伊逊河掌字水文站和不澄河边墙山水文站1990年7月14日至15日历时6小时流量过程 见图5 5 下游边界采用庙宫水库水位流量关系过程 见表5 1 验证计算中没有考虑排沙洞泄洪影响 假定大沽基面与国家高程基准面之间的换算关系为 假定大沽基面 16 868 1985国家高程基准 67 68 69 水流数学模型验证伊逊河掌水文字站和不澄河边墙山水文站1990年7月14日至15日历时6小时水位过程 验证结果见图5 8 验证计算流场见图5 9 70 71 从水位验证过程来看 计算值与实测值变化趋势一致 数值偏差不大 基本能保证模型预报计算精度 从不同时刻流场上看 当上游边界产生大流量时 公路桥附近收缩断面出现较大流速 下游大坝附近还没有达到溢流水位 当上游边界流量回落时 公路桥附近断面流速减弱 上游流量已传递到大坝附近 达到溢流水位 在坝前和库区弯道附近产生较大流速 基本符合水流运动规律 说明水流计算结果合理 水流模型可用于模型预报计算 从水位验证过程来看 计算值与实测值变化趋势一致 数值偏差不大 基本能保证模型预报计算精度 从不同时刻流场上看 当上游边界产生大流量时 公路桥附近收缩断面出现较大流速 下游大坝附近还没有达到溢流水位 当上游边界流量回落时 公路桥附近断面流速减弱 上游流量已传递到大坝附近 达到溢流水位 在坝前和库区弯道附近产生较大流速 基本符合水流运动规律 说明水流计算结果合理 水流模型可用于模型预报计算 72 2 泥沙模型验证泥沙数学模型验证采用计算模式验证的方法进行 验证1990年汛期前后平均滩面的变化 汛期上游来流过程见图5 6 上游边界含沙量过程见图5 10验证结果见图5 11 73 图5 10汛期上游边界含沙量过程 图5 111990年汛期前后水库平均滩面变化 74 5 4模型预报近年来随着水环境关注意识的提高和水资源高效利用思路的确立 如何在水库保水和水库排沙间选择优化的水库运行模式是需要进行研究和解决的重要问题 根据庙宫水库设计正常蓄水位778m 大沽高程 以下不加说明处均采用大沽高程 和库区死水位768m的要求 优选即能泄水排沙又能保水兴利的调水调沙方案是本次应用的主要目的 研究过程拟采用库区限蓄水位768m 770m 772m 774m 776m 778m六级水位标准 来水来沙过程采用1958 1959 1971 1975 1985 1989和1990的四个系列三种组合来模拟不同年份水沙过程 在18个计算系列中进行方案优选 下游排沙控制水位分别选取768m 770m 772m 774m 776m和778m 与三种水文类型年组合系列再组合形成18种计算方案 方案命名以水文类型年组合序号为先导字符 用横线连接排沙控制水位值 计算方案见表6 2 75 76 77 78 79 80 81 5 5排沙量和淤积量分析根据方案的计算结果 得到沿庙宫水库上游至下游各过流计算断面的床面平均变形量曲线 便于比较将同一水文系列的情况在同一图中给出 图5 33至图5 35 82 从图5 33至图5 35中可以看出 三种不同类型的来水来沙条件所塑造的河床变化形态相差较小 变化趋势完全一致 只是下游不同的排沙控制水位对河床冲淤趋势和冲淤量影响较大 故只给出水文条件下的最低排沙控制水位768以及最高排沙控制水位778情况下的库区含沙量 冲淤变化分布 流量沿程变化 见图5 36至图5 45 其他情况与此类似 应当注意的是 输沙含沙量与挟沙能力都是正值 而扬沙含沙量 冲淤量都存在正负 对于扬沙含沙量正值代表冲刷 负值代表淤积 冲淤量正值代表淤积 负值代表冲刷 83 84 85 86 87 88 89 90 91 92 93 94 95 从图5 33至5 35可以看出 不同的排沙控制水位 对于离坝距离大于5200m 铁路桥附近 的地方 影响很小 而影响范围 主要集中在大坝到铁路桥这一区段 变化趋势是 随着控制水位的不断增大 排沙洞到大坝这一段淤积增加 排沙洞附近冲刷能力加强 排沙洞上游淤积增加 总之就是 排沙控制水位高 离排沙洞较远的地方淤积就会加重 而排沙洞附近的冲刷反而会提高 进一步比较各方案之间在整体河道上的排沙效果 在表5 6中给出了庙宫水库在三种不同类型的来水来沙条件和不同排沙控制水位下各计算方案比较 13年水沙组合的来水量约为21470万 来沙量约为3174万 当下游排沙水位较低时 出库水量较大 出库沙量也较大 当下游排沙水位较高时 出库水量减小 出库沙量减少 可以看出 当下游排沙水位为水库死水位768m时排沙效果最好 在排沙水位为774m时仍有排沙效果 在排沙水位大于776m时水库出现淤积 当排沙水位取774m时即可起到排沙效果 又可蓄水兴利 为对各方案进行优选 设定水沙比和排沙比两个参数进行比较 其中 96 97 在表5 7中列出各方案的水沙比和排沙比 优选方案选择条件为 1水沙比尽可能小 2排沙比尽可能大 3排沙比不小于1 4排沙控制水位尽可能高 从表5 7中可以看出 1 774 2 774 3 774方案是在同类组合中水沙比最小 排沙比大于1 是符合上述四个条件的最佳选择 98 99 应用建立的数学模型和计算方法 通过Fortran数值计算软件编制计算程序 对庙宫水库库区内的水流和泥沙冲淤进行计算并与已有的实测资料对比 吻合效果较好 说明模型可以用于该水库的河道变形模拟 并且将该模型应用于不同水文年份 排沙控