（环境工程专业论文）太湖风生流及物质输运三维数值模拟.pdf

摘要 本文在对近年来同内外三维湖泊水流及物质输运数值模拟的有关文献进行系统分析 研究的基础上，进行了太湖风生流及物质输运三维数值模拟研究。论文从笛卡尔坐标系下 三维数学模型基本方程组出发，通过坐标变换，导出了垂向盯坐标系下的三维数学模型基 本方程组，给出了边界条件和初始条件；接着对基本方程进行了离散，采用过程分裂法求 解三维数学模型，并对紊动粘性系数进行了必要分析；论文对圆形水域和矩形水槽内的风 生环流两个理论算例进行了计算，得到了与解析解及前人计算较一致的结果，验证了模型 的精度；最后利用该紊流模型对太湖风生环流及物质输运进行了三维数值模拟研究。研究 表明在风应力的作用下，太湖水质同时存在水平和垂向变化，提出了太湖水质监测和评价 也要考虑污染指标的垂向变化。模型较好的反映了污染物质在太湖中的输移扩散规律，取 得了较满意的结果。该研究成果为今后建立完整的太湖水环境管理模型打下了良好的基 础。 关键词： 太湖三维 风生环流 物质输运数值模拟 a b s t r a c t b a s e do nt h ed e t a i l e da n da b r o a di n v e s t i g a t i o no ft h ei n l a n da n do v e r s e a sl i t e r a t u r e so f t h r e e - d i m e n s i o n a ll a k ef l o wa n dn l a s st r a n s p o r t a t i o nc o m p u t a t i o ni nr e c e n ty e a r s ，t h en u m e r i c a l s i m u l a t i o no f3 dw i n d d r i v e df l o wa n dm a s st r a n s p o r t a t i o ni nt a i h ul a k ei sb u i l t f i r s t ，t h e t u r b u l e n c ec o n t r o le q u a t i o n si n 仃c o o r d i n a t ei sc o n s t r u c t e df r o me q u a t i o n si nc a r t e s i a n c o o r d i n a t e ，a n d t h e b o u n d a r yc o n d i t i o n s a n di n i t i a lc o n d i t i o n sa l ea n a l y z e d t h e nt h e d i s c r e t i z a t i o ng r i da n dd i f f e r e n c es c h e m e sa r ep r e s e n t e d t h em o d e - s p l i t t i n gs o l u t i o nt e c h n i q u e i sa p p l i e di nt h i sm o d e l ，a n de d d yv i s c o s i t yc o e f f i c i e n t sa r ea n n y z e d b a s e do nt h e3 dn u m e r i c a l m o d e lt h ec o m p u t a t i o n so f w i n di n d u c e df l o wi nc i r c u l a rb a s i na n dac l o s e dr e c t a n g u l a rb a s i na r e t a k e n ，t h ec o m p u t e dr e s u l t sa r ei na g r e e m e n tw i t ha n a l y t i c a ls o l u t i o na n dt h er e s u l t so fl i t e r a t u r e s t op r o v et h i sm o d e lp r e c i s i o n f i n a l l y , b a s e do nt h em o d e l ，d i f f e r e n tw i n d - d r i v e df l o wa n dm a s s t r a n s p o r t a t i o ni nt a i h ul a k ea r et a k e n t a i h uw a t e rq u a l i t yi n h o r i z o n t a la n dv e r t i c a lv a r i e t y s i m u l t a n e i t yi sf o u n d ，v e r t i c a lv a r i e t yo fc o n t m n i n a t i o ni n d e xt ob ec o n s i d e r e di nt h el a k ef l u i d m a t t e rm o n i t o r sa n de v a l u a t i o ni sp u tf o r w a r d t h er u l e so ft r a n s p o r t a t i o na n dd i f f u s i o no f c o n t a m i n a t i o ni nt a i h ul a k ea r er e f l e c t e dp e r f e c t l y ，t h ec o m p u t e dr e s u l t sa r es a t i s f i e d t h es t u d y o ft h i sp a p e rm u s tb eu s e f u lt og u i d et h ew a t e re b v i r o n m e n tm a n a g e m e n tm o d e lo ft a i h ul a k ei n t h ef u t t 】r e k e y w o r d s ：t a i h ul a k e ，t h r e e d i m e n s i o n a l ， 、v i n d - d d v e df l o w ,m a s st r a n s p o r t a t i o n n u m e r i c a ls i m u | a t i o n 学位论文独创性声明： 本人所呈交的学位论文是我个人在导师指导下进行的研究工作 及取得的研究成果。尽我所知，除了文中特别加以标注和致谢的地方 外，论文中不包含其他人已经发表或撰写过的研究成果。与我一同工 作的同事对本研究所做的任何贡献均己在论文中作了明确的说明并 表示了谢意。如不实，本人负全部责任。 论文作者( 签名) ： 学位论文使用授权说明 年月日 河海大学、中国科学技术信息研究所、国家图书馆、中国学术期 刊( 光盘版) 电子杂志社有权保留本人所送交学位论文的复印件或电 子文档，可以采用影印、缩印或其他复制手段保存论文。本人电子文 档的内容和纸质论文的内容相一致。除在保密期内的保密论文外，允 许论文被查阅和借阅。论文全部或部分内容的公布( 包括刊登) 授权河 海大学研究生院办理。 论文作者( 签名) ：年月日 太湖风生流及物质输运三维数值模拟 1 1 研究背景及意义 第一章绪论 我国是世界上湖泊众多的国家之一，约有2 万多个湖泊，占全世界天然湖泊的1 1 0 ， 总面积约9 万多k 寸，其中1k w 以上的湖泊大约有2 7 6 0 个，湖泊总水域面积达7 万多k 种， 约占全国国土面积的0 9 5 。湖泊总贮水量达7 0 0 0 多亿m 3 ，其中约三分之一为浅水湖泊， 主要分布在东部沿海与长江中下游地区。由于这些湖泊地处平原，人口稠密，经济发达， 加上不适宣的人类活动方式，极大地改变了湖泊的天然生态系统的物质循环和能量活动， 以致于多数湖泊都面临着这样或那样的生态环境问题，如湖泊生态系统功能退化、蓝藻水 华频繁暴发、水质性缺水日趋严重、洪涝灾害加剧并造成巨大的经济损失等1 1 】。因此，以 太湖为典型，开展浅水湖泊生态环境的演化机理研究，对于满足国家在湖泊环境治理方面 的紧迫需求，有着非常重要的实际应用价值。 太湖是我国长江中下游地区著名的五大淡水湖之一。它位于长江三角洲的南翼坦荡的 太湖平原上，介于3 0 。5 5 4 0 ”3 1 。3 2 5 8 ”n 、1 1 9 。5 2 3 2 ”1 2 0 。3 6 1 0 ”e 之 问。地跨上海、江苏和浙江两省一市交界处，整个流域面积3 6 9 8 5k m 2 ，湖泊面积2 3 3 8k m 2 ， 湖泊蓄水量4 4 4 亿d ，居我国五大淡水湖的第三位。太湖南北长6 8 5 k m ，东西平均宽3 4 k m ， 最宽处5 6 k m ，平均水深1 9 m ，最大水深3 3 3 m ，是一个典型的浅水湖泊或积水洼地。湖底 地形十分平坦，平均坡度0 。0 1 9 6 6 ”，湖底平均高程为1 1 m ，地势由东向西倾斜，深 水区位于湖心略偏西，即大约在平台山以北水域。湖泊岛屿和暗礁很多，最大的岛屿是洞 庭西山。通常以洞挺东山的东胶咀为界，其西为广阔的太湖湖面，俗称西太湖，其东为东 北向的湖湾，俗称东太湖；气候上，太湖位于北亚热带，属季风湿润气候。四季分明，热 量充裕，雨水丰沛，光照充足。湖区年平均气温为1 5 3 1 6 ，年平均降水量为8 0 0 1 6 0 0 m m ，年蒸发量为8 0 0 1 0 0 0 r m n ，年日照时数为2 0 0 0 2 2 0 0 小时：太湖水源主要有三路： 南路来自苕溪水系，西路来自南溪水系，北路来自江南运河。其出水口多分布于湖的东部 和北部，主要有沙墩港、胥口、爪径口和太浦河等港渎下泄，再经吴凇江、黄浦江、望虞 ? 叭浏河等排入长江【2 一i 。 由此可见，太湖处r 我国经济最发达的睫江三角洲地区，是这一地区的重要淡水资源。 k 删以来，它已是上海、无锡、苏州等大、中城市的最重要的供水源地，且集供水、养殖、 旅游等多种功能于一体。因此。其水环境状况直接影响着这一地区的经济发展。但是，另 河海- 人学倾l 论文 一方面，由于近年来这一地区经济迅速发展，相应的环境保护和治理措施没有跟上，进入 太湖的j i 要河道和湖区的水质已r 益变坏特别是水体的寓营养化，已成为太湖水环境的 主要问题。按现行的国家地面水环境质量标准( g b 3 8 3 8 - - 2 0 0 2 ) ) ) ) ，太湖2 0 世纪6 0 年代 属i l i 类水体，8 0 年代初平均为i i i i i 类，8 0 年代末则全面进入i 类，局部i v 类和v 类，9 0 年代中期平均己达类，1 3 为v 类，到2 0 0 0 年太湖的2 3 3 8 k m 2 的湖面，i i i 类水质 占6 7 ( 1 5 6 7k m 2 ) ，i v 类水体占8 5 ( 1 9 8 9 7k m 2 ) ，v 类和劣于v 类水体占8 2 ( 1 9 1 6 k m ) 。太湖有8 3 5 的水域为富营养水平，仅1 6 5 为中营养水平，这主要在东太湖及湖 心地区，表现出太湖平均每l o 年左右，水质下降一个级别。2 0 世纪9 0 年代以后，蓝滦水 华每年暴发4 5 次，并逐渐向大太湖扩展。2 0 0 0 年监测结果显示，太湖湖心区域已经出 现严重的蓝藻水华。相应的，随着水体的富营养化，太湖的生物种群、组成也发生际变化， 蓝藻门中的微囊藻( 主要是铜绿微凛和水华微囊藻) 已成为藻类的优势种，浮游动物和底 栖动物大量减少，鱼类趋于小型化。种种现象都己表明，太湖的水环境变化正在加剧，其 生态系统自我调节功能和水体的自净能力正在逐渐衰退，太湖水质恶化造成城乡供水源的 污染和短缺，已成为该地区经济发展中新的制约因素，控制太湖富营养化和综合治理太湖 已是迫在眉睫。 要治理太湖，必须对太湖生态环境系统进行综合研究，才能达到规划、综合整治，这 就需要研究水、气、水中生物组成的自然生态系统，其中湖泊水动力学一直是湖泊研究的 重要方面，其主要动力过程湖流、增减水、波浪、定振波、深水湖泊的温跃层及内波 等对水体中各种物质、能量的传输及转化过程，水一气、水一土界面的物质交换存在着直 接和f h j 接的影响，是浅水湖泊水环境研究的基础。进行湖泊水动力过程的研究不但具有重 要的现实意义，而且具有重要的学科意义。业已展开对浅水湖泊水动力研究的结果显示水 动力作用是浅水湖泊生态环境演变极其重要的控制因子，对湖泊水域生态系统和环境有着 深远的影响。例如，富含营养盐的湖底沉积物在水动力( p c 、浪、湖流等) 作用下发生再 悬浮和迁移等，i j j 导致水体中悬浮质浓度增加。进而影响水色、透明度及初级生产力。太 湖藻类水华的暴发一般出现在大风作用过后的小风或静风、天气晴朗时段，这些现象说明 了水华的暴发与水动力过程有着某种联系。 国内外的湖泊专家早在5 0 年代就比较注意对湖流的观测。但由于受自然条件的影响， 获得湖泊波浪、湖流等要素长系列、高分辨率、在时间和空问上同步的资料难度很大，特 别是湖流受风场影响，稳定性差，其量级又低，观测十分困难。此外浅水湖泊水深一般都 太湖风生流及物质输运三维数值模扭 较浅，水体垂直方向动力混合作用十分强烈。所有这些，都加深了湖泊水动力学研究的难 度。进入7 0 年代，随着计算机技术及计算数学的发展，精细而全面地模拟水体运动过程 已能实现，再加上其耗资少，不受时间限制等优点，越来越受到重视。由加拿大内陆研究 中心、美国陆军团水道试验站、日本琵琶湖研究中心等研究单位的专家率先利用数学模型， 对北美的安大略湖( o n t a r i o ) 、伊利湖( e r i e ) 、日本的琵琶湖( b i w a ) 等世界著名湖泊的湖 流流忿进行数值模拟”l 。 本文主要针对当前太湖水质污染严重的问题，运用湖泊水动力学原理，建立湖泊三维 水流及物质输运数学模型，戍用该模型着力研究太湖湖流的特征、形成过程以及不同风速、 不同风向对太湖主要入湖河道排放口附近水域形成的污染物扩散范围的影响，力图找到一 些太湖湖流及入湖河道排殷口附近水域污染物扩散范围变化的基本规律以及具有使用价 值的分析结果，为今后建立完整的太湖水环境管理模型创造有利的条件。 1 2 湖泊水动力及物质输运数值模拟研究概况 湖泊水动力学一直是湖泊研究的重要方面，发展一个能够适用于湖泊水动力特性研究 的通用模型，从而推动湖泊环境学和水生态学研究的演化，是当前一个重要的热点课题。 日施太湖水动力学研究主要采用野外观测与数值模拟相结合方法，内容主要涉及风尘流、 吞吐流、波浪形成机理与分布特征，风浪、湖流作用下底泥再悬浮特征和对底泥营养盐释 放影i f 向，以及水动力作用下污染物与营养盐的输移扩散规律和风浪对藻类水华暴发的影响 等。 1 2 1 湖泊数学模型数值计算方法 湖泊水动力模型数值计算方法与其它水流的计算方法类似，目前较为流行的离散方法 有：有限差分法( f d m ) 、有限单元法( f e m ) 和有限体积法( f v m ) 、边界单元法( b e m ) 和分步法等。这几种方法各有优缺点及适用范围。 ( 1 ) 有限差分法( f d m ) 有限差分法适用与各种类型的偏微分方程，是目前应用最多最成功的一种方法，是以 尜恸级数展j f ：为工具将微分方程中的微分项以差商代替，从而化微分方程为易于求解的 代数方程。其特点是采用矩形网格概化计算域，各计算变量布簧在网格节点上。可分为显 ，法、隐式法、半显半隐法三类。在三维计算中，为了保证必要的精度并考虑到计算的经 济性，往往划空问的不同方向分别采用显式或隐式差分，其中交替方向隐式差分法( a d i ) 海人学坝上论义 p 1 是当前较为常用的一种差方法。有限差分法的优点是求解直观迅速，数学概念清晰，编 程及基础资料的整理比较简l 丫上，! 一用内存少，汁算的收敛性和稳定性到! 论比较成熟；其最 大的缺点是对复杂的地形适应弘l ：斧，网网格布置不灵活。难以概化不规则的边界。为了克 服其局限性，许多学者都致力于不规则边界处理问题的研究，如美国的t h o m p s o n j f 等 人提出的边界拟合坐标系( b o u n d a r yf i t t e dc o o r d i n a t es y s t e m ) 方法，l 谢利曼所提 出的任意网格差分法i ”，华祖林h 等也对边界拟合中的调节因子作过有益的探索。利用这些 方法在原则上可以把任意复杂的几何边界变成规则的几何边界求解，但计算区域的规则化 是以控制方程的复杂化为代价的，而且当边界存在尖角时，会出现局部奇异现象，使计算 不能收敛。 ( 2 ) 有限单元法( f e m ) 有限单元法也是种近似解法，是一种用积分求解微分方程的方法，从7 0 年代起开 始用于计算水力学中，主要是为了克服最初的有限差分法的规则网格划分不适应复杂边界 描述和刁j 能满足计算域内局部加密的要求等缺点而发展起来的。其原理是分单元刈觯迫 近，使微分方程空间积分的加权残差极小化。有限单元法具有网格划分灵活，拟合复杂边 界和地形容易等优点，程序通用性强，灵活性好，精度较高；其缺点是占用内存大，程序 编制比较复杂，对基础资料的整理比较繁琐，这些都限制了有限单元法的使用。另外随着 有限差分法在任意三角形和四边形的推广，有限单元法的优点已不占优势。常用的有限单 元法通常有r i t z 方法、g a e r k i n 方法以及迎风有限元法i ”。 ( 3 ) 有限体积法( f v m ) 有限体积法又称控制体积法，是7 0 年代由s p a l d i n g 和p a t a n k e r 等人提出和发展起 来的一种离散方法，最初应用于空气动力学领域。其基本思想是将计算区域划分成一系列 连续但互不重叠的控制体积，并使每个控制体积包围个网格点，以网格点上的因变量数 值为未知数，假定各变量在网格点间的变化规律，将待解的微分方程对每一个网格点积分， 得出一组离散方程，结合边界条件和初始条件求得数值解。有限体积法可以认为是一种结 合有限元方法改进的有限差分法，在网格点i 剞变量分如时，借鉴了有限元的思想，在离散 过程中应用的仍然是有限差分的方法。有眦体积法的优点是离散方程的解表示了变地( 如 质量、动量等) 对任意一组控制体积上的总体积分是守恒的，对整个计算区域自然也得到 守恒：有限体积法不仅能在精细的网格上获得较好的结果，而且在较粗网格的情况下，也 显示出准确的积分守恒。在网格布置时，往往采用交错网格，将不同的物理量布置在不同 的节点上此外，山于三维惬交网格的生成存在诸多困难，非诈交网格下的有限体积法i ”j 太湖风生流殷物质输运三维数值模拟 以及无结构网格下的控制体积法也得到一定的应用。因此，有限体积法在流体力学的数 值计算中得到了广泛的应用。 ( 4 ) 边界单元法( b e m ) 边界单元法早先称为边界积分方程法，在流体力学中称为有限奇点法或有限基本解 法。其实质是数学物理方程的格林函数法，据此给出解的积分表达式，再利用定解条件建 立边界积分方程。所得的边界积分方程通常不可能解析求解，因此利用有限元离散，将其 转化为边界节点一i ：未知量表示的代数方程组，求解该代数方程组，便求得边界节点所有的 未知值。由此可见，边界单元法可以看作格林函数法与有限单元法结合的产物。边界单元 法可以使求解的空间维数降低一维，从而使计算工作量和所需内存明显减少，对于三维水 流计算中自由表面的处理较为简单，其精度一般较有限元为高2 i 。这是边界单元法的最大 特点。边界单元法是计算椭圆性问题的有效方法，但由于需要控制方程的基本解，所以对 于复杂的问题，如完整的n s 方程尚未得到广泛的应用。 ( 5 ) 分步法 分布法，又称破开算子法，是将待解的偏微分方程分解为相对简单的若干部分，对于 不同部分采用不同的离散方法，建立不同的离散方程，用数值方法分别求解各部分方程， 从而得出偏微分方程的解m “】。分布法分为空间概念上的分布( 将多维问题分解成一系列 一维或低维问题来求解) 、物理概念上的分布( 将原始的包含若干个物理性质不同的复合 过程分解成若干性质单一的物理过程来求解) 和解析上的分布( 将原始解析方程进行分 解) - , f l , 类型。分布法求解微分方程包含两方面内容，一是如何分布，二是对各分布方程 采用何种离散方法求解。它是求解二维不流的种简便、易行而有效的数值计算技术，在 处理上具有相当的灵活性。原则上只要在分布之后，寻找出一个尽量正确、合理的数值离 散方式，就能取得良好的计算成果，分布法的主要不便之处在于处理边界条件。 在此基础上，湖泊水动力学模拟借鉴许多大气或海洋学中的先进方法，又发展了不规 则网格的有限差分数值模拟技术、有限分析法、边界分析法、嵌套网格技术等，使数值模 拟方法能更加贴近实际，渗透到了湖泊学研究的各个领域。 12 2 湖泊数学模型研究进展 国内外学者对湖泊水动力及物质输运模型的研究赋予了极大的热情，目前，常见的湖 洲数 = j ：f 型l 。j 分为三大类，即整层积分的：维模式、多层模式、三维模式。 ( 1 ) 二维模式 河海人学坝上论文 出于二维模式是将三维水流运动方程沿垂向积分后得到的，因此，二维模式又称积分 模式或全流模式，它适用于垂直混合强烈，不存在明显密度分层的浅水湖泊。它是由h a n s e n ( 1 9 5 6 ) 第一个提出的，用于计算浅水海域的水位变化过程及潮流，随后，h a n s e n 又用此 方法计算了康斯坦茨湖( l a k ec o n s t a n c e ) 环流模式( h u n t e r , 1 9 8 3 ) i ”j ；g a l l a g h e r l i g g e t t ，e ta 1 ( 1 9 7 3 ) 又建立了二维风生环流的数学模型【4 l ：在太湖水动力数值模拟研究方面， 吴坚( 1 9 8 6 ) 在国内首次设计了整层积分的二维浅水湖泊水动力学模式，模式采用了不规则 网格有限差分法，具有能确切地拟合湖岸、河口中、湖底地形及岛屿的功能，将浚模式应 用于太湖得到了太湖风生流、风涌增减水等有意义的结果，研究发现：太湖的l x t , q i 流极不 稳定，随风速变化而变化，并且稳定状态的湖流速较大，风涌增减水也十分明显，其幅度 取决于风力和风程f ；王谦谦( 1 9 8 7 ) 用二维差分模式模拟了定常作用下太湖风生流的各种 特征| | 7 j ；之后，王谦谦等于1 9 9 2 年建立了个嵌套网格的二维水动力模型刈太湖硐1 人汕w 口，风成流、风涌水作了数值模拟，结果表明：浚嵌套网格的模式对风浦水及流向有较好 的预报能力，并且能够节省大量的计算时| n i l ”i ：刘启峻( 1 9 9 3 ) 以守恒形式的二维水动力学 力脞组为丛砌，采刚。t7 a b a t a 的迎风柏限儿法模拟丫人湖，特别是太湖梅梁湾处赴各种风 况下的湖水运动过程，验证了太湖稳定流场分三种类型的正确性m 】；吴炳方等( 1 9 9 6 ) t i 二 维非恒定流模型分析了洞庭湖湖流运动规律和风力对湖流的影响作用i 1 i 。发展至今，二维 模式比较成熟，已广泛运用于实际工程的计算。二维模式的不足之处是不能得到速度的垂 向分布，对湖泊表面风应力及底部磨擦项常运用半经验关系。 ( 2 ) 多层模式 浚模式就是将水体沿水深方向( 垂向) 分成若干层，在每一层对方程积分，将三维问 题转化为若干个二维问题进行求解，因此该模式又称为分层模式。l e e n d e r t s e 等于1 9 7 3 年最早推出了这种模式，他假定各层中的变量在垂向上变化不大( 可以忽略) ，从而将三 维问题转化为多个二维问题，并在垂直方向采用固定分层法，即将计算水域分为固定的多 层，在每层中沿水深积分使之成为二维问题，并用a d i 格式进行数值离散l ：s i m o n 的大 面积湖泊环流多层模式( s i m o ne ta 1 ，1 9 7 3 ) ，在用于对安大略湖进行数值模拟时，将湖 体沿水深方向分为4 层，然后将各层分别处理为二维模式”j ：姜加虎( 1 9 9 1 ) 用一个二维分 层积分和一个二维整层积分的数值模拟方案，系统地探讨了抚仙湖、滇池环流形成机制， 分析了具有温度层结湖泊( 扰仙湖) 中内波以及上、下层流动的时空变化特征及其与湖区风 场的关系陋1 ：卢启苗采用分层二维的方法及平面集中质量有限元建立了海岸河口的三维潮 流数学模型1 ；白玉川在变换的垂向坐标下，采用分层二维法和算子分裂法，将三维控制 太湖风生流及物质输运三维数值模拟 方程分解为易于求解的物理意义单一的方程，用特征差分近似和o a l e r k i n 集中质量有限元 桐结合的方法对廉州湾潮流进行了三维模拟计算1 1 5j ；杨具瑞等将三维问题“二维化”的分 层迭代计算方法引入滇池湖泊水流计算中，对滇池进行湖流模拟，得到成功应用1 。 分层二维法简单实用，适用于实际分层流的计算，但该法在水深较大的水域计算比较 理想，在浅水部分的垂向分辨率很低，并且引入了人为交界面，交界面上必须给定摩阻系 数，摩阻系数的选用往往带有很大的经验性，且假定层与层之间是不可透水的，当流动存 在明显的上升或r 沉流动时，会导致计算的失稳。 ( 3 ) 三维模i 由于湖泊存雨直方向上也存在强烈的切变和混合作用，故许多的环境要素如溶解氧、 叶绿素、光强度及浊度等都存在有较大的垂直梯度，这就要求必须建立三维水动力模型来 模拟实际的湖泊水动力特征。d a v i e s ( 1 9 8 2 ，1 9 8 3 ) 开发了计算分层流的三维模式， s i m o n ( 1 9 8 5 ) 又，1 ：发了太湖风生环流的三维模式，澳大利亚水科学中心利用海洋水动力学 数值模式丌展了r 本琵琶湖的水动力学研究，k r e s t m tz i c ( 1 9 9 3 ) 建立了一个日本琵琶湖 三维斜压水动力学数值模式，详细模拟了湖泊环流结构及内波，但使用的风资料是湖中的 某一点的，不能代表湖中的实际风场f 4 i 。国内，焦春萌( 1 9 8 8 ) 提出了一个三维水动力学和 懋移物质输移的数学模型，模拟了太湖湖体中的动力过程和泥沙传输过程1 2 7 j ；逢勇等( 1 9 9 4 ) 采用浅水波方程数值模式，考虑太湖局地风场的影响，模拟了太湖风生流、风涌水的变化 情况，模拟结果与实测值符合较好；逢勇等( 1 9 9 6 ) 采用三维水动力学模式，研究了地形因 子对太湖环流的影响，指出洞庭西山和洞庭东山间狭长的过道区是产生太湖大范围环流的 霞要原因” ：王惠中等( 2 0 0 1 ) 建立了一个考虑垂直涡黏系数非均匀分布的太湖风生流准 三维数值模型，采用了交错网格技术布置各物理量和改进型的a d i 法，通过对稳定和不稳 定风场作用下太湖流场的模拟计算发现：垂直涡黏系数沿水深变化对水位的计算影响不 大，但刘流速的乖向分布和底床切应力的计算有较大的影响 2 。l 以上这些三维模型使人们对湖泊中的分层现象有了进一步的认识，但有时应用不是很 方便，例如在水深变化比较大的湖泊浅水区的层数可能减少。层面不同，模式的侧边界条 件位置不同，使计算程序的设计比较困难，使程序适用性受到一定限制，难以较好的反映 湖流的垂直结构。为了解决以上模型的缺点，一种垂直方向的压缩口坐标被从最初海洋研 究应用到湖泊动力学的数值研究。该方案1 9 8 3 年l = i = d a r a n c e 、t l u g h o s 等应用于欧洲北海m l ， 寞振兴川将盯坐标变换引入三维潮流计算中，并基于过程分裂法( 也称模态分裂法) 1 3 2 , 3 3 1 ， 即将三维流动中的快过程一表面重力长波与慢过程一缓行的内重力波分丌，并对每部分采 河海人学硕士论文 用适合其物理特性和数值行为的计算格式耦合求解得到了三维流场；1 9 9 4 年朱耀华把它应 用于我圈渤海、黄海、东海的海洋环境的研究1 3 4 1 ：粱瑞驹( 1 9 9 4 ) 应用了d r 模式训竹太湖 风生流i ；1 9 9 8 年胡维平应用了盯模式的完糇方程模拟9 7 1 1 号台风对太湖水位及流场的 影响( 胡维平等，1 9 9 8 a ) ，计算结果与6 点水位的观测值和单点流场过程符合较好；1 9 9 8 年胡维平等用此模型探讨了太湖各种不同风场对湖流的影响( 胡维平等，1 9 9 8 b ) ，比较了 整层平均流场与分层流场的差异，对浅水湖泊的二维模型作出了客观的评价 3 6 , 3 7 i 。 从理论上u l ：，三维模式自然更能反映实际的水流状态，尤其适用于对一些流域在垂直 方向非均匀分布湖泊的三维水流研究。但也存在一些不足，大规模的太湖数值模拟研究难 以开展，重复工作较多，造成了一定人力、物力、财力的浪费；计算花费时间较多，需要 比较完备的基础资料，而太湖基础性工作开展的较少，如对风拖曳系数和湖底磨擦系数的 实验确定等；再次，关于湖泊系统的实测资料较少，给数值模型的率定和验证带来了很大 的困难。 1 3 本文主要研究内容 目前，对于太湖风生流三维数值模拟研究前人已取得了一定的成果，湖流的计算已相 剥成熟，对风生流影响污染物质输运的数值模拟也有一定的研究成果，但研究主要集中在 二维模拟，运用三维模型研究太湖物质输运规律尚不多见。太湖风生流数值模拟的研究成 果表明，风应力引起的湖流各层差异很大，流速在垂向呈不均匀分布，从而影响污染物质 在垂向分布的不均匀性。因而，对太湖风生流及物质输运规律的研究更适合采用三维数值 模拟，本文的研究工作正是基于这一思路。 本文在总结前人工作的基础上，建立了湖泊三维风生流及物质输运数学模型，采用三 维水质模型模拟太湖污染物输移扩散问题。在水平方向采用笛卡尔坐标系，垂向采用 盯= ( x 。+ ) ( h + f ) 的坐标变换法，以利于提高浅水区的垂向分辨率，且垂向统一的分层 数也易于保证计算的稳定性，以提高计算的精度。本文着力研究太湖湖流的特征、形成过 程以及不同风速( 5 o m s 、2 5 m s ) 、不同风向( 只考虑典型风向东南风、两北风) 对太湖主 要入湖河道排放口附近水域形成的污染物扩散范围的影响，为今后建立完整的太湖水环境 管理模型创造有利的条件。 本文第一章对湖泊三维水流及物质输运数学模型的研究进展作出了分析和总结，提出 了本文的主要研究内容：第二章通过数学推导，导出了盯坐标系下紊流控制方程组，给出 太湖风生流发物质输运三维数值模拟 了相应的边界条件及初始条件，建立了a 坐标系下的湖泊三维水流及物质输运数学模型； 第三章采用有限差分法对仃坐标下的模型控制方程进行了数值离散，给出了内外模式方程 组的离散过程及差分格式，并对稳定性条件及紊动粘性系数进行了必要的分析；第四章利 用圆形水域风生流和矩形水槽风生环流两个典型算例对模型进行了验证；第五章对太湖风 生环流及物质输运规律进行了模拟计算，研究其湖流形态，并将其特点进行对比分析，同 时对太湖风生流的形成过程和太湖主要入湖河道排放口附近水域形成的污染带加以模拟 汁算；第六章对本文的工作进行了简要的总结，并指出了今后有待进一步深入研究的工作。 9 河海人学坝士论文 第二章湖泊三维水流及物质输运数学模型的建立 三维数学模型在维度上与实际水流一致，能尽可能真实的表现水流情况，近年来在河 l 二j 、湖泊水动力及物质输返刎：究中得到了j 泛的重视和研究。本文建立了适用于湖泊的三 维水流及物质输运数学模型，该模型具有以下的特点：在水平方向采用a r a k a w a “c ”i ” 网格，在垂直方向采用口坐标网格：采用盯垂向坐标变换后使整个水域具有相同的垂向分 层数，给数值计算带来便利；模型求解采用过程分裂法，将求解过程分为内模式和外模式， 即先求解二维外模式获得表面水位及垂向平均流速，接着求解内模式获得三维流场。下面 导出模型的基本方程组。 2 1笛卡尔坐标系下的三维数学模型基本方程 2 、1 ，1 常用假定 对天然水体直接求解三维n s 方程目前尚不可能，为此，人们提出了一些近似假定， 使实际天然水体三维问题具有求解的可能性，常用的假定有： ( 1 ) b o s s i n e s q 近似由于盐度、温度和物质浓度变化较小，除在动量方程的重力项 中必须考虑密度变化的影响外，其余各项一般可忽略不计，即在动量方程中，除重力项外， 视密度为常数。 ( 2 ) b o s s i n e s q 假定即b o s s i n e s q 在1 8 7 7 年提出的关于水流的紊动应力可类比于层 流粘性定力的假定。 ( 3 ) 静水压强假定天然湖泊大多数具有宽浅特征，垂向加速度远远小于重力加速 度，因而在垂向动量方程中忽略垂向加速度而近似假定为静水压强关系。目前水流运动的 三维数值模拟大多都采用该假定。 2 1 2 基本方程 基于以上假定，可以得到描述天然水域的三维水流运动的基本控制方程组。 连续性方程： 坐+ 旦+ 当：0 ( 2 1 ) o x l出，优1 动量方程： 0 太湖风生流及物质输运三维数值模拟 置枷毒+ v 毒+ w 差一一= i 去毒v 1 - + 毒( v t 塞vv + 杀。+ 毒乇，c z z ，百枷丙酉+ w 酉呻一瓦酉+ 刚喁j + 意t + 素乇， 旺2 1 r 一百i 产1 r 、， 卸却跏加。1 印 af 却1 aa 百棚百w 瓦+ w 酉+ 加一瓦毒+ 面l 坼两j + 意气：+ 索吃 旺3 j 0 2 一p o g ( 2 4 ) 融 物质输运方程： 筹+ “瓦o c + v 篆+ w 瓦a c = 毒p 篝 + 毒( 砧篆 + 毒( 厶鼍 + 氐 c ：s ， 紊流闭合模型；b - 程： 丝a t + “簧a x + v 堕a x 2 + w 簧= 旦a x 3 ( 生o k 十 薏a x + “+ v + w = 一| | - l 十r 。i 一1 鸭 i“j ，i + 毒 山誊 + 杀一盖 + 坼i 毒 2 + ( 毒 2i 一2 c z s ， 言c + w 毒c + v 杀c 盯，+ w 毒c 肼，= 三1 巨( 毒 2 + ( 毒) 2 + 毒卜毒c 肼， + 毒卜毒c 肼， + 毒 ( 尝+ 咋) 毒c 肼， _ 圭碥扩矿 c z ， 上述方程中：f 表示时间；薯，x 2 ，而分别表示笛卡尔坐标系中的x , y ，z 方向；“，v ，w 分别为 一，x z ，_ 坐标方向的流速分量；p 为压强：p o 为水体密度：c 为输运物质的浓度；s o 为源 项；k 为紊动能，为混合长；厂= 2 f i s i n # 表示柯氏力系数，q = 2 ，r 8 6 1 6 4 ：7 3 1 0 一r a d i s 表示地球的自转角速度；g 为重力加速度；嵋为紊动粘性系数；乃，砧分别为输运方程的 垂向和水平的物质扩散系数；a h 是七，k l 的水平扩散系数，本文取一= 0 2 v ；o k = 1 9 6 。 动量方程中，水平应力张量为： = 2 v h 丝0 x ih ：，叫萨h ( 豪+ 言) ，嘞毛妾 c ：s ，_ 1 一h z ，叫n 。h l 瓦+ 酉h z z - 2 蠢 乜” * u 婀人学删小论义 紊流闭合模型方程中，矿为壁面近似函数( w a l lp r o x i m i t yf u n c t i o n ) ，表达式为： 阳+ 巨阻去+ 去 2 ， 其中：盯= o 4 为卡f l 常数；e - e 2 = 1 2 盯2 ( 以岛。) ，e t ，e 2 为经验常数，一般取 e = 1 8 ，e 2 = 1 3 3 。 垂向紊动粘性系数叶和垂向扩散系数乃的表达式为： 坞= 最，七7 2 7 + ，五= s ，k v 2 z + 五 ( 2 ，1 0 ) 其中：和 表示背景混合系数，一般取v h = 以= 1 0 _ 6 肌2 s ：s o ，瓯为稳态函数( s t a b i i t y f u n c t i o n ) ，只。d 表示当n 2 = 0 ，即瓯= 0 时的稳态函数，由下列求解： 最= 蒜器，最= 器，q = 拿肚针盖考 。再j 而琵i 写研2 再面瓴2 i 几2 i l i 瓦j 眩 2 2盯坐标下数学模型基本方程 2 2 1 三维基本方程的推导 实际天然水体都具有自由水面，而自由水面的变动给水流的三维数值模拟带来了很大 的困难，本文采用盯坐标变换追踪自由水面1 3 0 l 。盯坐标变换的形式如下( 如图2 1 所示) ： 卜 h h 水面：z = 底面：z h 蚓2 1 承向仃坐标变换示意幽 ( i ，蜀，是，墨) = ( ：，矽( 盯) ) 计算区域 水面：0 = l ( 2 1 2 ) j i 一 概 一 0i1，llll l 太湖风生流及物质输运三维数值模拟 这里 仃：掣：型， 口【o ，i 】hhf + 。1 ( 2 1 3 ) 典- = i t 屯为笛书尔坐标系巾的垂向坐标；f 为自由水面相对于基准水位的波动高程，即水 位；7 , j 全水深，h 为平均水位到水底的距离，自由水面处z = 9 - ，水底面处z = 一h ：变换 后，墩f ( o ) = 0 ，( 1 ) = 1 ，当在自由水面时焉= l ，在水底部时墨= 0 ：雅可比变换为 j ：荽：j ，肌望 ( 2 1 4 ) 撕1ld o 2 2 1 1 笛卡尔坐标系与仃坐标系的转换关系 假设一坐标系中的时间为f ，- + ：h - 向坐标为量，五，墨，速度及物质浓度变量为 i ，i ，茹，e ，则笛卡尔坐标系与盯坐标系对应的坐标变换关系为4 1 i ： 文 兹 戚2 蜗 t = t 五= x 2 = x 2 毛= 巧( 口) ( 2 1 5 ) 为了推导速度变换公式，首先引入直角坐标对任意曲线坐标变换的全微分关系式： 蕊 威 幔 旅， 则关系式( 2 o o o0 c je e & a l获 氓& | 瓠。黩 瓠， a 曼| 8 ta 曼i f 瓠ia 曼t | 瓠：蕊t 瓠； 0 - 2 l | 8 t 蕊。| 钺l 蕊l | 瓠l赫：| 瓠3 舐，fa fa 曼3 瓠、a 曼； a x l蕊j a x 、 1 5 ) 对于仃坐标变换有如下的形式： 由上述坐标对应关系，可以推出： d f 如 如 如 出 豳 出2 凼 ( 2 1 6 ) ( 2 1 7 ) o o o l j 1llj 丝 仃一 o o o 一 鹭一 盯 ，一j 一 锄一 盯一 o l o - 望 口 一， j 塑! 坠苎兰堡主迨塞 于是得到速度转换关系 “= v = v 谛= 一筹；( 仃言- ( 1 一矿) 筹卜考( 仃亳仃) 剖+ 号 由于存在相等关系，以下的推导中将露，记做“，v ，从而导出垂向速度关系为 ( 2 1 8 ) ( 2 1 9 ) w = 国+ 盯熹十“( 嘻棚一盯，筹 + v 卜毒却仃，差) c z z 由关系式酉o f 2 瓦0 f + 嚣鲁和鲁= 毒( 口h 一一) = 盯等一筹 得到微分变换关系： 堡筹+ 蕊i a f l ( 孤a h o t一盯警 a f ，a i ，i 孤。a fj a f 舐l a f a 岛 o f 1 即f 搠阳1 酉十了酉l 酉可百j o fl a f l8 ha h 、 了1 _ = = - i = 一一盯一i 镌j 旎l 鹚。鹚 堡：堡堕：土堡 a ) c 38 曼3 融3 ja 曼1 ( 2 2 1 ) 6 h8 h8 8 8 h a h 百3 虿蔷。方百2 虿 帆a h0 0 e 嘲8 h i 2 鬲葛2 葛百2 面 2 2 ) 兰：旦，丝：堡塑：塑 钙一瓦 2 2 1 2 d 坐标系下连续方程的推导 由复合函数求偏导法则的关系式( 2 。2 2 ) ，笛卡尔坐标系下连续方程( 2 1 ) 中的各项 可以表达为： 詈2 詈与1 鹚9 ul ( 瓦。h 一仃筹 甄茂，蕊，i 蕊，。反j 4 ( 2 2 3 a ) 堕防一， + 峨万 堕 盯， 鸳一鸭 盯 一， 堕斫 孤一嘲 一 骘屯 ， m 、 一， 一 西舒武讲咄出口， | i = i i ，旗孑堕嚷万 太湖风生流及物质输运三维数值模拟 加加1a vf 舶o h1 酉2 瓦+ 了瓦l 瓦可酉j 跏1 挑 d x 3 jo x 3 由( 2 2 0 ) 式有 ( 2 2 3 b ) ( 2 2 3 c ) 篝2 专詈2 专毒l 西竹善卜誊吖) 筹j + v 卜薏- ( 1 一盯) 差 i瓦2 了酉2 了两i 抑竹蓄l 盯蔷一( 1 - 两j r 蠢一( 卜盯) 蠢j = 专去( 撕) + 吉喜+ 专卜瓦o h 一嚣 毒+ 号鼍+ 专( d 薏一薏 毒+ 号筹 ( 2 2 4 ) 将上述各式代入连续性方程( 2 1 ) ，可得： o u 1 锄f 劬o h l 西1 西f 劫8 - 1 1 ：+ 一? ij 一盯一i + 一十一一l 一一盯i 蕊l ，茂，l 蕊；撕喝，甄l 嘎鸥j + 专毒c 西，+ 吉等+ 专( 盯筹一嚣) 毒+ 号鬻+ 专( 盯瓦o h 一爱) 毒+ 青鼍= 。 ( 2 2 5 ) 整理后得到： 娶+ 孚+ 垫型+ 1 f 了o h + 。掣+ ，一3 1 - 11 ：o 酉+ 瓦+ 瞢+ 。i 虿州面w 瓦j = 0 j 百1 可o h + 吉( 詈+ “鼍 + 吉( 日薏+ v 薏 + 专掣= 。