（水文学及水资源专业论文）河道三维水流数值模拟.pdf

中文摘要 摘要 天然水流的实际运动呈现的是三维形态，因此要全面如实地反映水流的真实 特性，必须用三维的形式来模拟水流的运动状况，才能够更好地满足水利工程的 实际需要。本文对河道三维水流的模拟方式作了一定的研究，建立了能够应用于 实际的河道三维水流数学模型。 本文结合平面正交曲线网格坐标和垂向无因次坐标，建立了三维计算网格， 采用交错网格对各变量进行布置，结合冻结法对动边界进行处理，该方式可以方 便地处理露滩问题，并且可以较方便的实现程序编写。在所建立的三维计算网格 条件下对描述三维水流运动的方程组及其边界条件作了相应变换，对变换后的方 程算子剖开，第一分步方程组主要为流速求解，第二分步方程组主要为水位求解， 对破开后的方程组进行离散、差分，对流项采用“逆风”差分格式，扩散项用中 心差分，构造出沿主流方向的求解迭代面，在迭代面上采用垂向矩阵追赶法求解， 完成了河道三维水流数学模型的建立。在长江河道澄通河段进行模型验证，对模 型的参数( 糙率、垂向紊动粘性系数以及平面紊动粘性系数) 的选取进行了论证 说明，最终计算结果与实测资料拟合较好，证明该模型能够应用于河道三维水流 的模拟。 关键词：三维网格；算子分裂；迭代面；矩阵追赶法； 英文摘要 a b s t r a c t 1 1 1 ea c t u a lm o v e m e n to fn a m r a ln o wo ft h es h o wi s t h r e e d i m e n s i o n a ls h a p e ，t o a c c u r a t e l yr e n e c tt h eo v e r a l ln o wo fm et 】m ei d e n t i 坝m u s tb ei nt h ef o 咖o f n l r e e d i m e n s i o n a ls i m u l a t i o no ft l l em o v e m e n to fw a t e r ，i no r d e rt ob e 批rm e e tt l l e a c t u a ln e e d so fw a t e rc o n s e a i l c yp r 两e c t s i i lt l l i sp 印e r ，t i l r e e d i m e n s i o n a lf l o wo f t h er i v e rs i l n u l a t i o no fac e r t a i nw 批c a l lb eu s e dt oe s t a b l i s hm ea c t u a ln o w o ft h e r i v e rt l l r e e d i m e n s i o n a lm a t l l e m a t i c a lm o d e l n i sp 印e rc o m b i n ep l a i l eo n h o g o n a lc o o r d i n a t e 酊do fv e n i c a la n dn o n d i m e n s i o n a l c o o r d i n a t e s ，e s t a b l i s hat h r e e - d i m e n s i o n a lc o m p u t a b l eg r i d ，u s i n gt h es t a g g e r e d 鲥d l a y o u to ft l l ev a r i a b l e s ，t h ec o m b i n a t i o no ft h e 舭e z i n gp r o c e s so fm o v i n gb o u i l ( 1 a t l l i sm e t h o dc a i lb e e a s i l yd e a l t 谢t hb o t t o m l a n dp r o b l e m s ，a n dc a nb em o r e c o n v e n i e n tp r o c e d u r e sf o rm ep r 印a r a t i o no ft h er e a l i z a t i o n i nt e n n se s t a b l i s h e db y t h e 鹏e - d i m e n s i o n a lc o n l p u t a b l e 鲥du i l d e r 让怆c o n d i t i o no ft l 玳e d i m e n s i o n a l d e s c r i p t i o no f t 1 1 ec u 玎e n tc 锄p a i g i le q u a t i o n sa 1 1 db 0 1 1 n d a r ) ，c o n d i t i o n sc o 玎e s p o n d i n g 仃a n s f o 肌a t i o n ，t h e 觚l s f o 姗a t i o no ft h ee q u a t i o no p e r a t o rs p l i t ，t h ef i r s ts t e p - b y s t 印 e q u a t i o n sf o rm ef l o wr a t ec a l c u l a t em 旬o r ，a n dm es e c o n ds t 印_ b y s t 印e q u a t i o n sf o r t l l ew a t e r1 e v e lc a l c u l a t em 旬o r ，f o r c e dt h ee q u a t i o nf o rd i s c r e t e ，d i 虢r e n t i a l ，m e c o n v e c t i o ni t e mo ft ：h e ”u p w i n d ”d i 艉r c n t i a lf o m l a t ，廿l ed i f m s ei t e mo ft l l ec e n t e r d i f r e r e n t i a l ，a l o n gt h ed i r e c t i o no ft 1 1 em a i n s t r e 锄c o n s t n 】c t e di t e r a t i v es u r f a c e ，u s e m a t r i xc h a s e a r e rm e t h o dt o c a l c u l a t e ， c o m p l e t e d at h r e e d i m e n s i o n a lr i v e r m a t h e m a t i c a lm o d e lo fw a t e rn o w i i lt l l ey - 觚垂z e 慰v e rc h e n g t o n gr e a c ht og e t m o d e lv e r i f i c a t i o n ，d e m o n s 位旺em ep 姗e t e r so ft h em o d e l ( r o u g l m e s s ，v e n i c a l t u r b u l e n tv i s c o s i t ) ，锄dp l 锄et u r b u l e n tv i s c o s i t ) ，) s e l e c tw 1 1 i c hs h o wt 1 1 a tt l l ef i n a l r e s u h sf i tb e n e r 、 ，i t ht h em e a s u r e dd a t a ，n l a tt l l e l r e e d i m e n s i o n a lm o d e lc a nb e 印p l i e dt or i v e rn o w s i m u l a t i o n 1 ( q 唧o r d s ：t h r e e - d i m e n s i o n a l 鲥d ；批e t i co p e r a t o r ss p l i t ；i t e r a t i v es u 血c e ；m a t r i x c h a s e - a r e rm e t h o d ： 学位论文独创性声明： 本人所呈交的学位论文是我个人在导师指导下进行的研究工作 及取得的研究成果。尽我所知，除了文中特别加以标注和致谢的地方 外，论文中不包含其他人已经发表或撰写过的研究成果。与我一同工 作的同事对本研究所做的任何贡献均已在论文中作了明确的说明并 表示了谢意。如不实，本人负全部责任。 论文作者( 签名) ：。器二盔 毒，。年6 月侣日 学位论文使用授权说明 河海大学、中国科学技术信息研究所、国家图书馆、中国学术期 刊( 光盘版) 电子杂志社有权保留本人所送交学位论文的复印件或电 子文档，可以采用影印、缩印或其他复制手段保存论文。本人电子文 档的内容和纸质论文的内容相一致。除在保密期内的保密论文外，允 许论文被查阅和借阅。论文全部或部分内容的公布( 包括刊登) 授权河 海大学研究生院办理。 论文作者( 签名) ：二器盘 如8 年占月佃日 第一章绪论 1 1 问题的提出 第一章绪论 随着科技的进步和社会的发展，人们对各类流体的性质及运动规律的研究也 不断深入，对自然界和实际工程中的水体流动现象进行数值模拟是人们进行科学 研究的重要技术手段。计算流体力学是数值数学和计算机科学结合的产物，它以 电子计算机为工具，应用各种离散化的数学方法，对流体力学的各类问题进行数 值实验、计算机模拟和分析研究，以解决各种实际问题。 计算流体力学的基本特征是数值模拟和计算机实验，它从基本的物理定理出 发，在很大程度上替代了耗资巨大的流体动力学实验设备，该特点在科学研究和 工程技术中产生了巨大的影响。但由于水流运动的复杂性，在以往的科学研究中， 往往通过类比、假设及忽略某些次要因素等方法对真实水流进行模拟。按水流各 水力要素在空间上的分布特性来看，可以分为一维、二维和三维水流数学模型。 这三类数学模型各有自己的优缺点和适用范围【l 】：对于较平直、计算时段较长的 管道、隧洞、明渠等，一般只关心断面平均的水力要素，因此常采用一维水流数 学模型；对于河口、湖泊、宽浅河道等广阔水域地区，水平尺度远大于垂向尺度， 各水力参数( 如流速、水深等) 在垂直方向的变化远小于水平方向的变化，其流 态可用沿水深的平均流动量来表示，常采用平面二维水流数学模型；对于潮汐通 道、河口深渊段等窄深水域，相关参量的垂向变化远比水平横向变化为大，可以 采用垂向二维水流数学模型；二维水流数学模型基本上能够反映出水流的运动规 律，且二维模式简单，易于计算，得到了广泛应用。但是在海岸河口地区或受潮 汐影响明显的河道，水流的水平尺度与垂向尺度较为接近，当研究水流运动、河 床演变、水质污染扩散等问题时，二维数学模型就不能很好的满足实际要求，应 该采用三维数学模型，其模拟效果好，但控制条件和编程计算复杂，计算成本高。 总体而言，2 0 世纪7 0 年代以一维水流计算为主，进入8 0 年代后，除继续应 用推广一维水流计算外，大多已采用二维水流数学模型，并根据实际工程需要， 配合以泥沙、温度、盐度和污染物等物质输移模型。9 0 年代以来，随着工程项目 的需要，三维水流模型的研究和应用日趋广泛，并取得了一些颇有价值的研究成 河海大学硕士学位论文 果。但由于其复杂性，目前往往采用简单的模型来代替复杂的数学模型，即用在 一个或两个方向上的积分平均，将三维的问题降为多个二维或一维问题以简化三 维流动结构。 但是，水流的实际运动呈现的是三维形态，要全面如实地反映水流的真实特 性，必须用三维的形式来模拟水流的运动。与二维水流数学模型相比，三维模式 虽然计算量大，但是能很好的模拟水流的空间结构，更精确地反映出水流的真实 运动规律，能够更好地满足水利工程的实际需要，而且，随着计算机的普及与进 步、计算机性能的提高，为三维水流数学模型的建立提供了强有力的技术支持， 因此，如何更好地建立三维水流数学模型代表了今后的水流数值模拟技术的发展 方向，有很好的研究意义和现实意义。 1 2 三维水流数值模拟的研究现状及存在问题 从1 9 3 3 年，1 1 l o m 采用手算的方法得到了不可压缩粘性液体绕圆柱流动的差分 数值解开始，三维粘性流动的数值计算方法在流体运动的研究中得到了越来越多 的关注【2 】，从建立数学模型、给定定解条件、生成计算网格、选择数值计算方法 等诸方面作了大量而细致的研究工作，并且一些成果已经在实际工程中得以应 用，从而促进了三维流动数值模拟技术得发展，尤其在潮汐河口中的水流运动呈 现的三维特征非常明显【3 】，具有很高的研究价值。因此，近几十年来，很多学者 对我国河口及其邻近海域的三维水流运动提出了基本数值模式和计算方法。 1 9 8 5 年，赵士清【4 】提出了一个三维数学模型，能够应用于长江口外海滨潮流 的计算，该模型在垂向上水体被划分为若干层，从而将一个三维问题简化为若干 个二维问题，各层之间通过内摩擦阻力来连接，基本方程采用显式有限差分来近 似，基本单元为正方形，对空间变量的偏导数采用中心差分，对于时间变量的偏 导数采用显式差分，对于水深较大的水域，用该模式进行三维潮流计算，结果比 较理想。 1 9 8 9 年，韩国其等【5 l 以静水压强假定下的不可压缩三维非恒定雷诺方程和 茁一g 双方程湍流模式为基本方程，模拟充分混合的潮汐河口中三维流速分布和水 位随时问的变化过程，引入无尺度垂向坐标，使计算网格能拟合自由表面和床面， 对微分方程利用破开算子法求解，通过使用欠松弛技术提高了计算的稳定性，对 2 第一章绪论 长江口潮汐水流的计算表明，模拟能真实地再现潮汐河口中的实际流动。周思平， 张书农1 6 j 针对明渠中三维底孔喷射流问题，提出了一个新的求解s 方程的有 限差分法( h h s i m p l e 方法) ，用该方法求解三维水流问题，可以保持 s i m p l e 方法的优点，即可将动量方程中的速度与压力解耦，求解离散化动量方 程时可用隐式迭代格式，而无需作“刚盖假定 或“静水压力假设 ，因而可以 很好地处理自由表面。 1 9 9 0 年，宋元平、胡方西1 7 j 根据长江口外盐度分层现象明显的特征，建立了 一个三维数学模型，此模型是建立在一系列二维数学模型的基础上，通过将水体 在垂向上分成若干层，在每层水体上对流速进行平均处理，来计算长江口外海滨 的潮流和盐度扩散变化情况，侧重研究了盐度分布在垂向上的变化。 1 9 9 2 年，徐贵泉等【8 】在周思平的博士论文h h s i m p l e 方法的基础 上，对其补充方程二的离散格式，作了进一步改进，使得改进后的h h s i m p l e 计算模式可适用于一般的潮汐河口的往复流动和更复杂的三维含自由水面的流 动，从而扩充了该模式的应用范围，并针对长江口北槽情况，也采用r g 紊流模 型来闭合三维非恒定雷诺方程，建立了三维数学模型，进行了数值计算。1 9 9 2 年 曹德明，朱耀华等【9 】和1 9 9 3 年李身铎，顾思美应用由方闰洪和朱耀华发展的一 种仃坐标三维非线性模式，在水平方向上以较小尺度的差分网格覆盖计算区，垂 直方向上给予均匀的分层，较好的模拟了杭州湾的潮波运动。 1 9 9 6 年，刘子龙，王船海等【5 5 】针对当前具有自由表面的三维非恒定水流数 值模拟所存在的问题进行分析，用平面正交曲线坐标变换结合垂向无量纲变换， 将整个计算区域变换成固定的规则区域，然后选择合适的算子分裂，根据每一分 步算子的特性，构造各自合适的求解方法，来解决大型方程组的求解问题，紊流 模型采用简单的零方程模型，以此编制了一个较通用的二、三维连接计算程序， 对长江口南支吴淞口附近的实测资料进行验证，计算值与实测值非常吻合，并且 该程序具有露滩处理功能，计算精度较高，工作量小等特点。 1 9 9 7 年，华祖林，卞华【1 1 】针对开敞性水域三维水流的特点，抛开了垂直粘性 系数的假定，而以静压假定建立垂向水深无量纲水流运动三维计算模式，引入了 r 一占紊流模型求解垂向粘性系数 ，通过连云港海域实例计算，较好地模拟了水 流在垂向上的变化分布。王晓建，张廷芳【1 2 】由三维浅水方程出发，利用有限体积 3 河海大学硕士学位论文 法的积分模式并引入g r e e n g a u s s 公式，将控制体单元的体积分转化为面积分， 采用一种任意六面体单元，垂直涡动粘性系数定义为水深的函数，使之能适应复 杂海域和地形，又可获得不同水深的计算数据，用该方法对辽东湾海区流场进行 了数值模拟，获得了理想的结果。 2 0 0 0 年，华祖林【l3 】采用拟合曲线网格系统模拟天然河道边界的基础上，建立 了符合弯道( 含急弯河道) 水流特性的贴体坐标系下三维数学模型，并针对弯曲 河段可能引起的回流及环流，引入了双方程茁一g 紊流模型精细模拟，经长江某急 弯河段的实例计算，取得了较为理想的效果。李缇来，窦希萍【1 4 】在平面上用正交 曲线网格、垂向上用等间距分层、变量交错布置，建立了长江口边界拟合坐标的 三维潮流数学模型，较好地反映出长江口的河口形态和潮流运动。 2 0 0 4 年，李艳红，周华君【l5 】将平面正交曲线坐标系和垂向盯伸缩坐标系相结 合，采用模式分裂技术将控制方程分解为外、内模式联合求解，首先由外模式计 算出水位梯度和水深平均流速，再由内模式向外模式提供流速值，以计算外模式 中的底部切应力，如此反复循环计算，最终使方程收敛，得出水位及流速分布的 正确结果，对方程的离散采用有限差分法，建立了适合弯曲河流的三维水流数值 模型，对嘉陵江广元城区河段进行模型验证，成功地模拟出了水面横比降和弯道 二次流，模拟结果与实测数据吻合较好。 2 0 0 5 年，吴修广，沈永明等【1 6 j 对c a n e s i a l l 坐标系下地凡蝌s 方程进行三维 善一刁一f 坐标变换，建立了非正交三维曲线坐标下弯曲河流的标准k g 湍流模 型，通过非正交曲线坐标下平面二维河流数学模型对计算区域进行二维计算，求 得全场水位，用这一水位作为计算区域的表明，沿这一曲面作为“刚盖 来近似 自由水面，与天然河床高程划分三维数值网格，在河床阻力的模拟中，采用壁面 函数法来近似，用控制体积法离散方程，该模型通过具有实验数据的实验室连续 弯曲水槽进行验证，模拟的流速值与实验数据吻合良好，并将模型应用于天然连 续弯曲河流的流场计算，显示该模型具有模拟天然弯曲河流的能力。 2 0 0 6 年，赖锡军，曲卓杰等【1 。7 】建立了新型的基于特征高分辨率数值算法的三 维非结构网格浅水动力模型，该模型采用有限体积法离散s i 舯a 坐标下的三维浅 水方程，运用r o e 黎曼近似解评估水平界面通量，模型网格拟合边界能力强，可 根据需要局部加密，格式数值性能优良，具有守恒性，单调迎风性，高数值分辨 4 第一章绪论 率等特性，同时，应用干湿判别法处理动边界，以适应浅滩地形漫露过程模拟的 需要，封闭水池内部风生环流，干河床上溃坝过程和长江口实际潮流场的模拟从 不同侧面展示了该模型的特点，结果表明它能够准确地预测水流的三维流动结 构，而且计算简单高效，具有良好的数值稳定性。张明亮，沈永明等【1 8 】采用非正 交三维曲线坐标下的彭一占双方程湍流数学模型，对弯曲不规则地形的渠道和漫滩 流动进行数值模拟，该模型在水平方向采用非正交曲线坐标，在垂直方向采用等 分网格的全坐标变换，对河道底部表面和两岸壁面使用壁面函数法处理，在划分 网格时，对漫滩和主槽分别单独划分曲线网格，再合成统一的曲线网格进行计算， 对顺直漫滩带有弯曲主槽的漫滩流动和大型水槽进行了详实的三维计算，结果表 明该模型对三维漫滩水流流动有一定的应用价值。 2 0 0 7 年，董壮，龚政等【1 9 】针对仃坐标变换的不足，结合绝对分层的方法，引 入垂向坐标变换系数，通过严密的数学推导，建立了改进的三维水流传质数学模 型，并给出了改进后的控制方程及相应的边界条件、初始条件，采用有限差分法 对模型进行了离散，针对殴江河口局部海域进行了潮流、盐度过程的验证计算， 取得了较好的成果。 综上所述，近三十年来，三维水流数学模型的研究取得了显著进展，已经逐 渐成为了研究水流实际运动状态的强有力工具。但同时，上面介绍的多种三维数 值模拟方法的应用仍然具有局限性，仍需要进一步的深入研究。如目前应用较多 的盯坐标变换，盯坐标变换的优势在于对自由水面的处理相对简单，但同时也存 在一定不足，即变换后的控制方程形式更为复杂，计算结果需加以转换，计算难 度增大，还可能引入新的误差，某些情况下容易发散，在地形变化复杂处误差较 大等等【2 0 】。 1 3 本文研究方法与内容 鉴于目前三维水流研究方法的多样性，查阅了近几十年来的一些三维水流数 值模拟方法，分析了它们的特点，建立了三维水流数学模型，即在平面上采用正 交曲线拟合坐标，在垂向上采用无量纲变换，并对控制方程和边界条件进行相应 的变化，对动边界采用冻结法处理，使得对动边界的计算简化成固定边界的计算； 对变换后的控制方程用算子分裂法求解，构造了沿主流方向的求解迭代面，在迭 5 河海大学硕士学位论文 代面上采用垂向矩阵追赶求解，完成了模型的建立，并使用该模型在长江河道澄 通河段进行了模型验证，得到较好的效果。本文的主要工作如下： 第一章绪论。叙述了选题的目的和意义，指出了研究三维水流数学模型的必 要性和重要性，并对近几十年来三维水流数值模拟的研究现状进行了综述，提出 了本文的研究方法与内容； 第二章河道三维水流数值模拟的研究方法。指出了描述水流运动的控制方 程，并概述了简化方程常采用的假定、对自由水面的处理、常用的网格生成方法， 主要对控制方程的一般离散方法作了简单介绍； 第三章河道三维水流计算模型的建立。主要介绍了描述三维水流运动的基本 方程组及其边界条件，着重讲述三维计算网格的生成和对三维水流方程的离散差 分及求解。采用平面正交曲线网格和垂向无量纲变换，进行节点布置和垂向网格 划分，以此构造了三维计算网格，并介绍了对动边界的处理方法；对差分后的方 程组中的对流项采用“逆风”差分格式，扩散项用中心差分，采用算子分裂法求 解，构造了沿主流方向的求解迭代面，在迭代面上采用垂向矩阵追赶法求解，并 设计了程序的计算流程，完成了河道三维水流数学模型的建立，并经过恒定流测 试，来证明本模型的正确性、可行性； 第四章模型的验证应用第三章所建立的河道三维水流数学模型，经过反复 率定确定了糙率、垂向紊动粘性系数等影响模型计算结果的参数，然后对长江河 道澄通河段进行了模型验证，结合实测潮位站点资料及实测流速测点资料，对模 型的计算结果进行分析，证明该模型能较好地反映水流的三维特性，能模拟河道 三维水流运动； 第五章总结与展望论述了主要研究工作，存在问题及进一步研究的建议； 6 第二章河道三维水流数值模拟的研究方法 第二章河道三维水流数值模拟的研究方法 2 1 描述水流运动的控制方程 不可压缩粘性流体的质量守恒方程，即连续性方程为： 丝+ 堕+ 堕：o ( 2 1 ) 瓠 动 a z 式中：，甜，分别为x ，y ，z 方向的速度分量。 假设为所研究水体的通用变量，则按普遍的守恒性原理满足【2 1 】： 詈+ 昙c 力+ 参c 力+ 妄c 力一昙c 巨等一号c 髟爹一鲁c e 警，= s ( 2 2 ) 式中：t 为时间，“，“：分别为x ，y ，z 方向的速度分量，矽可以代表速 度分量畋，或水质指标的浓度等；巨，髟，疋为综合扩散系数，s 则为 相应于通用变量的源项。 方程( 2 1 ) 与方程( 2 2 ) 构成了封闭的方程组，从理论上来讲，结合适当的 初始条件和边界条件，求解上述方程组就可以得到对水流运动的真实描述。但由 于该方程组是复杂的非线性偏微分方程组，直接进行数值求解有一定的难度，于 是，人们根据问题本身的物理特性，提出了一系列的假定与近似，以使实际问题 得到简化。 2 2 简化方程常采用的假定 2 2 1 常用假定 1 b o u s s i n e s q 近似 自然界中的水体，其温度、盐度和物质浓度随空间和时间变化，这必将引起 水体密度的变化。但在实际问题中，水体密度的变化相对较小，且对质量力产生 7 河海大学硕士学位论文 的影响也很小，因此，在动量方程中，除了重力项考虑密度变化外，其他各项一 般可忽略不计。即在动量方程中，除重力项外，视密度为常数。这一近似已被广 泛应用于大体积水体密度流和分层流的计算中。 2 b o u s s i n e s q 假定 由于在紊流时均运动方程组中包含了较难处理的雷诺应力张量，b o u s s i n e s q 在1 8 8 7 年提出关于将水流的紊流应力( 雷诺应力) 类比于层流的粘性应力，与 时均速度梯度成正比的假定，即用层流粘性应力的形式对紊动应力进行参数化。 3 静水压强假定 天然水流大多具有宽浅特性，垂向加速度远远小于重力加速度，因此在垂向 动量方程中忽略垂向加速度而近似假定为静水压强关系。但是对于复杂突变边界 附近的水体，水流状态急剧弯曲、升降或回流，此时若采用静水压强假定会带来 较大的误差。 2 2 2 自由水面的处理 自然界中广泛存在具有自由水面的水体，如何在三维水流数学模型中追踪模 拟自由表面是研究的重点之一，也是难点之一，通常用以下几种处理办法【1 】： 1 刚盖假定 假设在自由面上存在一个固定的刚性盖，限制自由面位置随时间的变化，从 而使沿垂向的网格剖分及计算易于进行。实践证明，该假定对计算时段内水位变 幅不大时，既能反映流动的主要特征，又能使计算得到合理的简化。但是，如果 研究自由表面水流的局部变化或自由水面变动较为剧烈的情况下，则刚盖假定不 能取得令人满意的结果。目前，对刚盖假定进行改进是一个主要的研究方向。 2 标记网格法( m a c ) 该方法采用一组随流体运动的标记粒子，这些粒子不含质量，本身不直接参 与计算，用只有坐标位置的标记粒子来描述自由表面所在位置圆，在某一时刻， 含该粒子的区域为流体域，不含该粒子的区域为非流体域，在两域的交接处即为 自由液面。该方法可以较好地模拟含自由表面的不可压缩流体运动，但收敛性较 差。 3 体积率法( v o f ) 8 第二章河道三维水流数值模拟的研究方法 v o f 法改变m a c 法中对全部流场用标记点进行标记的做法，只对自由表面 进行跟踪。基本思想是：定义一个体积率函数f ，表示计算单元内某种液体体积 占计算单元体积的相对比率。当计算单元完全充满液体时，f = l ；当计算单元为 空时，f = o ；当f 在0 和1 之间时表示该计算单元被液体部分充满，区域内存在 自由液面。计算时首先确定哪些单元包含有自由水面，然后根据f 函数的梯度来 确定自由水面的法线方向。计算出各单元的f 数值及其梯度之后，就可以确定各 单元中自由水面的近似位置。v o f 法相比m a c 法可大量节省内存和计算时间。 5 拉格朗日不可压缩法( l i n c ) l i n c 法采用拉格朗日网格，不采用标记点( m a r k e r ) ，也不采用质点( ( p 础i c l e ) 。 网格节点随流体运动，每个单元中的流体永远保留在本单元中，流体的边界随着 单元边界运动。l i n c 法可用于计算带有自由表面的粘性不可压缩流体的非恒定 运动，较易考虑表面张力效应，适合于求解流体变形不大的问题，否则网格将发 生交叉、扭曲等现象而使计算无法进行。 2 3 常用网格生成方法 网格生成是数值求解中的一个重要组成部分，也就是对空间上连续的计算区 域进行剖分，把它分割成许多个子区域，并确定各个区域中的节点。网格生成的 方法很多，陶文铨【2 3 1 在其论著中作了总结( 见图2 1 ) 。网格生成主要分为三大类： 结构化网格、块结构化网格和非结构化网格。 9 河海大学硕士学位论文 网 格 生 成 方 法 结构化网格 正交曲线坐标系中的常规网格 适体坐标法 对角直角坐标法 块结构化网格 二二二 保角变化法 代数法 微分方程法 撇格臣： 吖： 图( 2 1 ) 网格生成方法分类 结构化网格单元一般为四边形或六面体，每一个节点及控制体积的信息必须 存储，但该节点的邻点关系可以根据网格编号的规律自动得出，而不必专门存储 这一类信息。其网格节点有序排列，网格中的任一节点与周围节点的关系固定不 变。结构化网格具有结构简单，节点信息存储方便存储量小的优点，而且以结构 化网格为基础己发展了许多成熟高效的离散格式和算法。例如传统的有限差分法 和有限体积法在求解流体动力学问题时均采用直角坐标系下的方形网格或矩形 l o 法 法 法 法 方 程 程 程 化 法方 方 方 范 法 值型 型 型 规 界 法 插圆 曲 物 界 边 面 限椭 双 抛 边 双 多 无 ril fj、l 第二章河道三维水流数值模拟的研究方法 网格。这对于简单而规则的流动区域是合理的，可以保证网格的光滑性和正交性。 但实际上许多流动现象发生在不规则区域，直角坐标系下的方形网格或矩形网格 划分就给方程求解带来诸多不便。目前，在复杂边界流动计算中，适体坐标网格 是应用比较广泛的结构化网格。 采用块结构化网格时，计算区域需要分解成两个或两个以上的子区域。在每 个子区域中均用一种形式的结构化网格，各区域之间可以互不重叠，也可以互相 重叠。在重叠区域，信息的传递是研究计算方法的一个重要方面。块结构化网格 具有结构化网格的特点，同时不要求一条网格线贯穿在整个计算区域中，给处理 不规则计算区域带来不少方便。 非结构化网格一般为三角形或四面体，节点的排列灵活，节点的编号与命名 无一定规则，甚至是完全随意的，而且每个节点的邻点个数也不是固定不变的， 每个节点与相邻节点的连接方式也不尽相同，非结构化网格表现出一种不规则、 无固定结构的特点。与结构化网格相比较，非结构化网格的优点是网格结构灵活， 对复杂区域有很好的适应性，网格生成的自动化程度高，容易实现自适应算法， 缺点是缺乏与之适应的高精度的离散格式，数据结构复杂，占用内存多；另外由于 数据结构的复杂性，求解代数方程的很多经典有效算法，如线迭代方法、a d i 方 法、块修正方法在非结构化网格中难以应用。 2 4 控制方程的离散方法 对水流控制方程的离散主要有【2 4 1 ：时间项的离散、空间项的离散( 包括对流 项、扩散项) 。计算流体力学里常讲的离散方法主要针对对流项的离散。随着计 算机性能的提高及计算方法的发展，离散方法朝着高精度、高分辨率、高效的方 向发展，归纳起来主要有一下几种： 1 有限差分法( f d m ) 有限差分法用有限差分方程逼近微分方程，求解差分方程组得到离散点上物 理量的差分近似解【2 5 1 。如果差分方程组不耦合求解，即各方程可独立求解，称为 显格式，反之，若需联立求解，称为隐格式。差分格式的分类标准很多，例如按 逼近精度的阶分为一阶、二阶，以至更高阶；按格式的性质分为中心及逆风( 或偏 心) 格式两大类。根据离散格式的不同，有限差分法应用较为广泛的是以下两种求 河海火学硕士学位论文 解方法： ( 1 ) 直接差分法 直接差分法是将微分方程中的微商直接用差商来代替，它可分为显式差分格 式和隐式差分格式。显式差分格式可以直接求出未知节点上的值，不需要联合求 解，然而显式差分格式收敛性和稳定性较差，时间步长不能过大。隐式差分格式 克服了显式差分格式稳定性和收敛性差、时间步长受限制、精度低的缺点，在数 值计算中有广泛的应用。 ( 2 ) 特征线法 特征线法基于沿特征线成立的特征方程( 又称相容关系) ，而不是利用普通 空间坐标中的原始方程进行流体流动的计算。特征方程反映了双曲问题中信息沿 特征线传播的性质，因而算法符合水流的物理机制2 6 1 。从本质上说，特征线法是 一种合理的逆风格式，具有优良的精度。不足之处是，特征方程常为非散度形式 ( 非守恒形式) ，用差分法离散特征方程时会带来较大的守恒误差、节点扰动大， 当处理多维问题时，过于烦琐。 总之，有限差分法建立在经典的数学逼近理论的基础上，简单且易为人们接 受，处理效率较高，是发展较为成熟的方法，在计算流体力学方面获得了广泛的 应用。差分方法的计算精度较其它一些方法稍低，从理论上讲可以通过提高差分 的阶数来提高精度，但又会带来其它不便，因此，目前主要通过改进差分格式来 解决这一问题。差分方法的另一个弱点是，运用差分方法的矩形网格不能较好地 适应复杂的边界。针对这一问题，近年来提出了多种解决办法，如运用梯形、三 角形等单元拟合边界，但应用较为普遍的是适体坐标法。 2 有限单元法( f e m ) 有限单元法将计算域划分为若干个单元，在各单元内部选定适当的节点，将 微分方程中的变量改写成由各变量或其导数的节点值与所选用的插值函数组成 的线性表达式，应用变分法或加权余量法，生成有限元离散方程，将各单元的离 散方程合成为总体方程组，结合初值、边值条件得到近似解。有限单元法主要优 点在于处理不规则边界或域中局部点函数值变化显著的问题【2 7 】。 由于控制方程的变分形式不存在，计算流体力学中应用较为广泛的是加权余 量法，选择不同的权函数就可以得到不同的方法，如子区域法、最小二乘法、 1 2 第二章河道三维水流数值模拟的研究方法 g a l e r 飚n 方法等。在传统的有限元方法基础上还有多种改进方法，如质量集中有 限元、迎风有限元等。隐式有限元方法的缺点是对非恒定流计算时，每一个时间 步长要求解一个大型线性方程组，耗机时多，计算所需存储量很大：同时，常用 的g a l e r k j n 有限单元法在数学上适于求解椭圆型方程的边值问题，不适于求解以 对流为主的输运问题：还有，g a l e 凼n 有限单元法的性能类似中心差分格式，不适 于计算间断。 3 控制体积法( f v m ) 在差分法的基础上，发展了控制体积法。控制体积法又称有限体积法，其基 本原理是将计算域划分为一系列不重复的单元( 称为控制体积) ，以网格点上的因 变量为未知数，并假设其在网格点之间的分布规律，将控制方程对每个控制体积 积分，合成一组离散方程，结合初边界条件求得数值解。从积分区域的选取方法 来看，控制体积法属于加权余量法中的子区域法；从未知解的近似方法来看，控 制体积法属于采用局部近似的离散方法。简言之，子区域法加离散，就是控制体 积法的基本方法【2 引。 4 边界单元法( b e m ) 边界单元法是以边界积分方程为基础，用有限元思想将边界积分方程离散 化，导出近似的代数方程组。其基本做法是将近似解的剩余用满足微分方程的奇 异函数( 格林函数) 为权做内积，分步积分直到微分算子全部转移到权函数上，推 出与微分问题等价的边界积分方程，然后在边界上剖分单元，离散后得到代数方 程组。因此边界元法是一种以格林函数为权的特殊形式的加权剩余法。 边界元法有许多优点，首先它以边界积分方程为基础，只需在边界上求解， 所以能使问题的维数降低一维。如原为三维空间的微分问题可化为边界面上的二 维积分问题，而原二维问题可以降为一维问题。它只需将边界离散化，而不需将 整个区域离散化，减少了代数方程系数矩阵的阶次，减少了计算工作量及所需的 计算机容量，但是边界元法得到的是稠密的方程系数矩阵，与有限元法得到的是 稀疏矩阵相比，求解带来诸多不便；同时需要己知所求解偏微分方程的格林函数 基本解，虽然对不少偏微分方程的格林函数基本解已经找出，但对s 方程这 样的非线性偏微分方程，至今未找到其基本解，目前只能采取s 方程中的非 线性项看作扩散方程的源项并通过迭代的方式来求解低雷诺数的流动。 河海大学硕士学位论文 5 有限分析法( f a m ) 有限分析法【2 9 】【3 0 】【3 q 是陈景仁教授于1 9 8 0 年提出的一种新的数值方法。其基 本思想是在网格单元上求解局部线性化的非线性偏微分方程，由微分方程的解析 解构造离散代数方程。它的离散格式保留了原微分方程的特点，能够根据对流的 方向和大小自动调整格式的系数，避免了在差分近似中的各种数值效应。有限分 析法具有自动迎风、精度高、稳定性好等优点，因而得到较为广泛的应用。但由 于有限分析系数中含有无穷级数，给实际计算及理论分析都带来了一些困难；同时 计算工作量较大，对计算区域几何形状的适应性也相对较差。 另外常用的控制方程的离散方法是分步法，分步法又称算子分裂法，是由苏 联学者y a n e n o k 刚和m a r c h u k 【3 加等提出的，其基本思想是将原来比较复杂或高维 的数学、物理问题的求解分解成若干比较简单问题的连续求解过程。分步法是求 解多维问题中构造比较经济的数值格式的一种有效方法，在工程界得到普遍应用 【3 0 】f 3 4 】【3 5 】【3 6 1 。分步法有三种类型，空间概念上的分步、物理概念上的分步和解析 上的分步。在空间概念上分步的破开算子法是目前流行的一种数值计算方法，其 基本思想是将一个微分算子分裂成多个简单算子的线性组合，在各个分步步长上 分别求这些简单算子的解，从而得到整个步长的数值解。 p e a c e m a l l ，r a c h f b r d 和d o u 翊a s 于1 9 5 5 年提出了a d i 法( 显隐交替法) ，该 方法把时间步长t 分为两个半步长，前半个时间步长化在x 方向用隐格式， 在y 方向用显格式，后半个时间步长化在y 方向用隐格式，在x 方向用显格 式，这样可以把较大的系数矩阵转化为两个三对角矩阵，用追赶法求解。a d i 法 不像全隐式那样每一步都要解全计算域上的大型方程组，在计算成本上比较经 济；同时，它在计算过程中一般比较稳定，因此，该方法在提出之后很快得到了 推广和应用。在此基础上，许多学者又进行了一系列改进，发展了一系列新方法： 例如1 9 6 7 年l e e n d e n s e 【3 7 1 提出了显式与隐式的混合模式，采用a d i 求解二维明 渠水流，为了克服非线性项的影响导致格式呈现不稳定的问题，对原格式进行改 进，采用分步全隐式，形成了l e e n d e r t s e 法；还有水位、水深、流速等变量均布 置在同一个网格节点上的m a d i 法( m o v i n ga d im e t l l o d ) ，该方法既吸收了原有 a d i 算法的优点，又有较高的稳定性、收敛性和精度。 1 4 第二章河道三维水流数值模拟的研究方法 2 5 本章小结 本章指出了描述河道三维水流运动的控制方程，并概述了简化方程常采用的 假定，尤其是对自由水面的多种处理方式作了介绍，然后介绍了常用的网格生成 方法，主要对控制方程的一般离散方法作了介绍。对有限差分法、有限单元法、 控制体积法等作了简单介绍后，对算子分裂法以及a d i 法作了较为细致的介绍。 1 5 第三章河道三维水流计算模型的建立 第三章河道三维水流计算模型的建立 3 1 三维水流的基本方程组及其边界条件 3 1 1 基本方程组 在满足第二章所述的假定条件下， 程： 连续方程： 锄加却，、 + + = 0 舐砂 钯 得到笛卡尔坐标系下的三维水流动力学方 ( 3 1 ) 明重力程： 票+ “罢+ v 罢+ w 罢一声+ g 等：昙p 为+ 吴 罢) + 昙( t 当 ( 3 2 ) + “+ v + w 一砂+ 2 = i 占) + ( g ) + l 占l j = z ) 国反 却瑟。苏叙、衙7 却、却7a z 一瑟7 安+ 甜祟+ v 安+ w 妻+ + g 票：昙p 尝) + 昙( 占安) + 昙( t 罢) ( 3 3 ) + 甜+ v + w + 础+ _ g 一= ( f ) + i 占) + ( s ) l 3 3 ) 街缸 却 a z 。 。却舐、叙7 却、却7 a ! z a ! z 7 式中：x 、y 、z 为空间笛卡尔坐标，u 、v 、w 分别为三个坐标方向上的流速， 压为水位，占、乞分别为水体的平面及垂向紊动粘性系数，厂为柯氏力参数，取 7 r = 2 q s i n 矽，q 为地球旋转角速度，够为地理纬度，g 为重力加速度。 3 1 2 边界条件 3 1 2 1 边界条件概述 初始条件和边界条件是确定微分方程解的必不可少的条件，初始条件产生的 误差在计算稳定的条件下可以逐步消除，而边界条件的误差会直接影响整个计算 的精度。河道水流计算中，可作为边界条件的有潮位、流速和流量。在有选择的 情况下，通常应选水位作为边界条件较为适宜，因为水位沿横断面的变化往往较 小，即所谓的横比降较小；而且边界应选在地形比较平缓和流态为缓流的地方， 1 7 最好选在有流速、潮位资料的断面上，并且断面应与流线保持正交。 3 1 2 2 边界条件确定 1 自由表。面边界条件 一胁胁竺 ( 3 4 ) 比2 百枷s i + k 万 “1 互= ，= 风应力( 对j 求和) ( 3 5 ) 2 两岸边界 互= m 乃= 岸边摩阻力( 对j 求和) ( 3 6 ) 可丙= 0 ( 不穿透条件) ( 3 7 ) 3 底部边界 嘞= 鲁饥鲁+ 屹誓 。8 ) z = ，= 底部摩阻力( 对j 求和) ( 3 9 ) 4 水与水交界处的边界条件 水与水交界处的边界条件分两种情况： 矽丙为已知，水位乏已知； ( 3 1 0 ) 孑丙为已知，流量q 己知，同时假定该断面水位无横比降； ( 3 1 1 ) 上面各式中：压为水位，乙为河底高程，秒表示流速矢量，= ( l ，2 ，3 ) 为边界面的法向矢量，z ( 待1 ，2 ，3 ) 表示边界处的切应力3 8 1 ，叱表示自由水面流 速，w d 表示河底髓摩阻力采用斌昭磐可瓤拶西，鼽c 与n 分别表示谢才系数与糙率。 1 8 第三章河道三维水流计算模型的建立 3 2 三维计算网格的生成 3 2 1 概述 流体数值计算是在离散的计算网格上进行的，因此，采用的计算网格是否适 宜，对计算的成败起关键作用，好的计算网格应该能完整、准确地反映边界的情 况，疏密适当，既能较好的控制计算量，使计算有较好的收敛性、经济性，又能 使计算结果误差