（流体力学专业论文）近岸海域沉积物示踪模型及其应用.pdf

中文摘要 近年来，随着人们环境保护意识的提高，近岸海域的污染排放有所减少，但 是对于历史上沉积于底泥中的重金属、p a h 和d d t 等污染物，由于其难分解、 残留时间长，很容易随着底泥的启动造成水体的二次污染，破坏生态环境，甚至 危害人类健康。因此研究近岸海域沉积物的迁移扩散规律对于对近海污染的防治 有非常重要的意义。 本文首先建立了近岸海域的水动力学模型( 包含波浪模型和潮流模型) 、污 染物迁移扩散模型和泥沙模型，并在此基础上，提出了一种理想示踪剂的概念， 建立了沉积物示踪模型。对于波浪模型，首先用经典的实验结果进行了验证，然 后根据渤海湾的实测资料对其进行了校正和验证。对于潮流模型，用渤海湾海域 实测数据对模型进行了验证，模拟结果与实测结果符合良好。对于泥沙模型，分 别计算了潮流单独作用下以及潮流和波浪共同作用下的泥沙浓度分布；根据曹妃 甸附近的大潮平均含沙量对纯潮流作用情况做了验证，结果吻合较好；对于风浪 天气下的渤海湾含沙量情况，用2 0 0 0 年9 月6 日和2 0 0 1 年3 月8 日的渤海湾遥 感含沙量分布图进行了验证，尽管在具体数值上有所差别，但模拟的泥沙运动趋 势是一致的。各项验证结果表明本文建立的水动力学、对流扩散模型和泥沙数学 模型是可靠的，可以支持沉积物示踪模型的运行。 用验证后的水动力模型计算了渤海湾的余流场，通过余流场的特性对沉积物 的可能迁移路径进行了分析。然后模拟了释放于渤海湾海域的保守物质的迁移扩 散行为，初步得出了保守物质的迁移扩散路径。最后用本文建立的沉积物示踪模 型模拟了沉积物的迁移行为，得出了沉积物的迁移扩散的一些基本结论。结果表 明：保守物质和沉积物的长期迁移路径与余流的特性是一致的；由于扩散的效应， 保守物质和沉积物除随余流迁移外，还不断的向外海扩散；随水深的加大，悬浮 泥沙将逐渐下沉，并且不易再次起悬，因此，与保守物质相比，沉积物更难迁移 到外海。 关键词：沉积物：示踪剂；数学模型；渤海湾；波浪；潮流 a bs t r a c t i nr e c e n ty e a r s ，w i t ht h ei n c r e a s eo ft h ea w a r e n e s so fe n v i r o n m e n t a lp r o t e c t i o n ， t h ep o l l u t a n td i s c h a r g e di n t oc o a s t a lw a t e r sh a sd e c r e a s e d h o w e v e r , s o m ep o l l u t a n t s ， s u c ha sh e a v ym e t a l s ，p a h ，d d ta n ds oo nc a ne x i s ti nt h ec o a s t a lz o n e sf o rav e r y l o n gt i m e ，d u et ot h e ya r ev e r yd i f f i c u l tt ob ed e c o m p o s e d t h e s ep o l l u t a n t sa r ee a s l y t ob ea d s o r b e dt os e d i e m n t sa n dd e p o s i ti nt h es e ab e d u n d e rs o m ec o n d i t i o n s t h e s e p o l l u t a n t sc a ne n t e ri n t ow a t e rb o d ya g a i nw i t ht h er e s u s p e n s i o no ft h eb e ds e d i m e n t ， a n dd a m a g et h ee c o l o g i c a le n v i r o n m e n t t h u s ，i ti sv e r yi m p o r t a n tf o rt h ec o a s t a l w a t e re n v i r o n m e n t a lp r o t e c t i o nt os t u d yt h eb e h a v i o u ro ft h es e d i m e n t s t h eh y d r o d y n a m i cm o d e l s ，w h i c hi n c l u d e sa2 dw a v em o d e la n da2 dt i d e m o d e lw e r es e tu pa tf i r s t t h e nt h ea d v e c t i v e - d i f f u s i o nm o d e la n dt h es e d i m e n t t r a n s p o r tm o d e lw e r ee s t a b l i s h e db a s e do nt h eh y d r o d y n a m i cm o d e l s f i n a l l y ，t h e c o n c e p to fi d e a lt r a c e rw a sp r o p o s e da n dt h es e d i m e n tt r a c e rm o d e lw a sd e v e l o p e d t h eh y d r o d y n a m i cm o d e l sa n dt h es e d i m e n tm o d e lw e r ev e r i f i e da g a i n s tt h e e x p e r i m e n t a ld a t a ，t h ef i e l ds u r v e yd a t aa n dt h er e m o t es e n s i n gd a t a t h er e s u l t ss h o w t h a tt h em o d e l sd e v e l o p e di n t h i st h e s i sc a nb ea p p l i e dt os i m u l a t et h et i d a lf l o w ， w a t e rw a v ea n ds e d i m e n tt r a n s p o r ti nb o h a ib a y t h er e s i d u a lf l o wf i e l dw a sc a l c u l a t e du s i n gt h eh y d r o d y n a m i cm o d e la n dt h e p o s s i b l et r a n s p o r tp a t ho ft h es e d i m e n tw a se s t i m a t e db ya n a l y z i n gt h er e s i d u a lf l o w f i e l d t h e nt h eb e h a v i o ro ft h ep a s s i v ec o n s e r v a t i v em a t e r i a li ss i m u l a t e du s i n gt h e a d v e c t i v e d i f f u s i o nm o d e l f i n a l l y , t h et r a n s p o r to ft h es e d i m e n tw a ss i m u l a t e du s i n g t h es e d i m e n tt r a c e rm o d e ld e v e l o p e di nt h i sr e s e a r c h t h en u m e r i c a lr e s u l t ss h o wt h a t a ) t h el o n gt i m et r e n do fs e d i m e n tt r a n s p o r ti sd o m i n a t e db yt h er e s i d u l af l o w b ) t h e m o v e m e n to ft h es e d i m e n ti sd i f f e r e n tf r o mt h a to ft h ep a s s i v ec o n s e r v a t i v em a t e r i a l d u et ot h es e d i m e n tw i l ld e p o s i tt ot h eb e du n d e rc e r t a i nh y d r o d y n a m i cc o n d i t i o n s ， a n dc ) i ti sd i f f i c u l tf o rt h es e d i m e n t ( a n dt h ep o l l u t a n t sa d s o r b e dt ot h es e d i m e n t ) t o t r a n s p o r ta n dd i f f u s ei n t of a rf i e l d k e yw o r d s ：s e d i m e n t ；t r a c e r ；m a t h e m a t i c a lm o d e l ；b o h a ib a y ；t i d e ；w a v e 独创性声明 本人声明所呈交的学位论文是本人在导师指导下进行的研究工作和取得的 研究成果，除了文中特别加以标注和致谢之处外，论文中不包含其他人已经发表 或撰写过的研究成果，也不包含为获得叁鲞盘堂或其他教育机构的学位或证 书而使用过的材料。与我一同工作的同志对本研究所做的任何贡献均已在论文中 作了明确的说明并表示了谢意。 学位论文作者签名：兰8 鹃年l 签字日期： 勘。c i 年6 月5 - 日 学位论文版权使用授权书 本学位论文作者完全了解墨鲞盘鲎有关保留、使用学位论文的规定。 特授权鑫鲞盘堂可以将学位论文的全部或部分内容编入有关数据库进行检 索，并采用影印、缩印或扫描等复制手段保存、汇编以供查阅和借阅。同意学校 向国家有关部门或机构送交论文的复印件和磁盘。 ( 保密的学位论文在解密后适用本授权说明) 学位论文作者签名：兰胭穹每三 签字日期：2 。u c j 年6 月f 日 导师签名：么心、签 签字吼刀年多月f 日 第一章绪论 1 1研究背景及意义 第一章绪论 我国是一个海域辽阔，海上利益广泛的国家。海洋对我们来说，是寄予很大 希望的、具有战略意义的资源接替空间。海湾、近岸等海域是外海和陆地的交汇 地带，是人类活动频繁的场所，并与沿海地区居民的生产和生活息息相关。海洋 渔业、养殖、海上交通、运输、海洋石油开采和海滨旅游等，都对沿海经济的发 展与人民生活的进步有着重要的作用。但是近年来，随着经济的发展和人类对海 洋开发的不断深入，近海海域的生态环境负担越来越大，海域环境质量明显下降， 生态环境日趋恶化，并对生物资源和人体健康产生有害影响。近海环境保护与污 染治理问题成为全球关注的焦点问题之一。 由于固体( 包括底泥和悬浮物) 颗粒的表面吸附效应，使海底沉积物成为很 多排入海水中的污染物的主要藏身之处，它能从天然海水中富集污染物，可以说 是海水中污染物的净化剂。但天然海水的环境条件变化会改变与沉积物结合的污 染物的可迁移性，导致污染物重新释放出来，再次进入生态环境，造成二次污染， 从而危害人类健康，破坏生态环境。对于一些持久性的污染物，如重金属、p a h 和d d t 等，由于沉积物的储蓄和释放引起的二次污染效应尤为显著。近年来， 随着环保意识日益深入人心和环保措施的逐步到位，近岸海域的污染排放大大减 少。随着外源输入性污染物的比重逐渐降低，底泥释放的污染物的比重逐渐增大， 因此，研究近岸海域沉积物的迁移扩散，对污染的防治有非常重要的意义。 渤海作为我国唯一的内海，南北跨度5 5 0 公里，东西宽3 4 6 公里，海岸线长3 7 0 0 公里，面积约7 7 万平方公里，渤海海峡宽仅1 0 6 公里，典型的腹大口小。平均水 深1 8 米，最大深度8 6 米，是一个深入大陆内地的近封闭的浅海。根据国家海洋局 公布的资料显示，n 2 0 0 0 年底渤海海域4 0 的区域已经受到严重污染，整个海域 呈现富营养化状态。其中，严重污染海域主要分布在人口密度大、工业区集中的 大中城市沿海地区。渤海湾是渤海三大海湾之一，是一个典型的半封闭海湾。由 于邻近京津等工业城市，陆源排污量大，海洋污染严重。本文选取渤海湾天津海 岸段作为研究对象，就是因为渤海湾海岸带在环渤海各地区中相对最脆弱、退化 最严重，生境修复需求突出，且该段海域历史上曾受到比较严重的污染( 包括重 金属、p a h 和d d t 等持久性污染物) ，沉积物中污染物的含量较高，底泥释放造 第一章绪论 成的二次污染对水体总污染的贡献率可能较大。 1 2国内外研究现状 河口海岸沉积物运动是风，浪，流共同作用的结果，建立余流的沉积物迁移 转化及其示踪模型，必须考虑到以上各项动力因素。 1 2 1波浪 通常，近海波浪是指风波、涌浪和混合浪，三者皆属于短重力波。对于波浪 场数值模拟的研究经历了由规则波到不规则波，由线性波到非线性波的过程，人 们提出了许多关于波浪传播与迁移的计算方法，由于每种方法的物理假定不同， 其适用范围也不同。近海常用的波浪模型有b o u s s i n e s q 类方程模型、缓坡方程模 型和s w a n 模式。其中前两者主要用于涌浪的模拟，不适宜于大面积海域风浪 场的模拟；而基于波作用量守恒的s w a n 模式在近年来广泛用于近海风浪场的 模拟。 8 0 年代后期，s w a m pg r o u p 将适用于深海的相平均谱波浪模型进行改进， 提出了适用于近岸及沿岸风浪预报的动谱平衡方程。之后，t o l m a n 雎1 提出了基于 动谱平衡方程的w a v e w a t c h 第三代波浪模型；青岛海洋大学1 等单位也进行 了海浪数值预报研究；r i s ，h o l t h u i j s e n 和b o o i jh 儿副总结有关来波浪能量输入、损 耗及转换的研究成果，提出并发展了适用于海岸、湖泊及河口地区的s w a n 模型 ( 第三代浅水波浪数值预报模型) ，全面合理地考虑了波浪浅化、折射、底摩擦、 破碎、白浪、风能输入及波浪非线性效应。b o o i j 等以线性波浪理论及现场观测 资料对s w a n 模型进行了验证，结果表明，第三代浅水谱波浪模型能比较准确合 理地模拟复杂潮流、地形、风场环境下的波浪场。 目前广泛使用的s w a n 模型( 含源代码) 由荷兰d e l f t 大学土木工程系开发并 维护。从第一个公开发布的版本s w a n 3 0 5 l 到目前的4 0 5 l 版本，进行了逐步的 改进和扩充，性能不断提高，功能也不断增强。它采用了基于e u l e r j 匠似的波作 用量谱平衡方程，采用线性随机表面重力波理论( 包含流的作用) 哺3 。该模型考虑 了较多的物理过程，包含了当前海浪预报研究的最新成果。 王殿志口引，张娜【1 和郝宇驰5 】等用s w a n 模型对渤海湾的风浪场进行了数值 模拟。其研究结果表明，在对渤海这种大范围浅水波浪场进行模拟时，s w a n 模 型是合理的选择。 第一章绪论 1 2 2潮流 潮流是近海和河口地区最主要水动力条件之一，是这些地区物质运动的最基 本形式，其他物质诸如泥沙、盐分、各类污染物质及热量的输移，均伴随着潮流 而运动。 潮流数值模拟始于上世纪6 0 年代，7 0 年代末便有大量研究成果问世。从控制 方程来说，潮流场数值模型可以分为二维模型和三维模型。三维模型因其计算精 度高、可模拟潮流与泥沙的垂向分布，越来越受到人们的关注并已逐步形成了一 批较为成熟的计算模式，被应用于区域较大的天然流场的计算中口引。而在河口、 海湾等水域，由于空间水平尺度远大于垂向尺度，在底坡比较平缓的情况下，可 以对垂向的运动进行概化，通过对三维潮流运动进行垂向积分得到平面二维潮流 模型。二维潮流模型因其计算量小，在其适用范围内应用广泛。李世森嘞3 等应用 变网格有限元方法，采用四边形等参单元离散横向积分的n a v i e r s t o k e s 方程，建 立了简单有效的河口潮流垂向二维数学模型。江春波等哺”将分步有限元法引入浅 水方程的求解，并采用了一些增加计算稳定性的因素，开发了适用于模拟具有复 杂边界浅水流动的计算模式。孙健舱2 1 、孙涛乜副和何磊幢4 1 等用水深积分的平面二维 模型对渤海和渤海湾的潮流场进行了模拟，表明对于渤海湾，二维模型可以获得 满意的水位和流速计算结果。 1 2 3沉积物 对于海湾沿岸来说，最主要的沉积物就是泥沙，因此研究海湾沿岸的沉积物 情况主要就是研究泥沙的相关情况。目前，研究水沙问题的方法主要有两种：物 理模型和数学模型。这两种方法都是以现场资料作为试验的基础，但各有其适用 范围和优缺点。就物理模型而言，主要存在着比尺效应，投资大，周期长，可移 植性差，很难完全适应多因素、大范围、多方案的工程规划问题。而物理模型的 缺点正好是数学模型的优点。近几十年来，随着计算机与计算技术的发展，数学 模型被越来越广泛地应用于解决近岸工程实际问题。 泥沙数学模型的研究是从一维开始的，经过数十年的研究和应用，一维泥沙 数学模型已经比较成熟。在国内，诸裕良瞳引等建立了一维河网非恒定流及悬沙数 学模型，提出了节点悬沙控制法，并成功用于珠江三角洲流域河网水流及悬沙的 数值计算。林秉南瞳础等通过分析泥沙沉积速度与含沙量之间的关系以及不同时段 半潮平均挟沙力，建立了钱塘江河口的一维悬沙模型。二维泥沙模型可以对泥沙 的对流、扩散、冲刷、淤积等动力过程以及底部泥沙固结机制进行模拟，二维泥 沙模型一般可分为平面二维和垂向二维模型。平面二维泥沙数学模型建立在垂向 第一章绪论 平均的基础上，能够模拟出区域泥沙场的平面分布，而垂向二维能够反映泥沙剖 面的分布以及泥沙底边界的某些动力过程。m e a d 瞳引同时使用垂向二维和平面二 维的泥沙模型研究河口内挖槽区的淤积问题，并对其做了比较。在我国，利用二 维数学模型来研究泥沙问题发展也很迅速，如张士奇等曙踟利用一、二维衔接的数 学模型对黄河口冲淤进行了研究，在河道中用一维模式，在河口用二维模式。叶 锦培等心蜘对涨落潮使用不同挟沙力公式，建立了二维泥沙模型并应用于珠江口磨 刀门潮流输沙的计算。曹文洪等啪1 开发建立了适合黄河河口近岸应用的平面二维 动边界非恒定水流泥沙数学模型。严以新等呤门建立了三维水域全隐水沙数学模型 并成功应用于珠江三角洲航道整治工程。窦国仁等2 1 建立了波浪和潮流作用下包 括悬沙和底沙在内的全沙数学模型。随着计算机技术的迅速发展，多种数值模拟 软件也逐渐出现，如丹麦d h i 开发的m i k e 2 1 ，m i k e 3 和荷兰d e l f t h y d r a u l i c 开发 的d e l f t 3 d 等均可较好的模拟沉积物运动。 对于渤海湾沉积物的分布及迁移情况，以及渤海的环流、潮余流对沉积物的 影响作用，赵宝仁等m 1 ，王辉、郑连远等7 儿3 踟在9 0 年代通过监测数据或数值模拟 结果得出了一些结果，但各家对余流的结论不尽相同，赵宝仁等根据实测海流资 料数据得到的渤海中部及辽东湾内的平均环流是顺时针的，而王辉、郑连远用数 值方法得出的环流则是呈逆时针运动的。赵宝仁还根据实测资料得出了渤海的平 均环流与渤海沉积物的分布和输运相协调的结论。 1 2 4示踪技术 所谓示踪技术，就是从外面加入与所需跟踪的对象完全共同运行的示踪 物，用以追踪所研究对象的运动或变化的一种方法。所谓示踪剂，就是人们为 观察、研究和测量某物质在指定过程中的行为或性质而加入的一种标记物。常 见的示踪剂有同位素示踪剂、酶标示踪剂、荧光标记示踪剂、自旋标记示踪剂 等。 示踪技术的应用领域非常广泛，在水文地质等工作中的用途都很大。人工 示踪剂是人类有意识引入的非天然成分的示踪剂，是一种非常有效的研究手 段；这种方法可以将运移参数量化，同时也可以很好地刻画水环境特征。另外， 采用示踪剂可以直接进行原位测试，选择物理化学性质合适的示踪剂，还可以 研究某些非常特殊的作用过程。在许多情况下，示踪剂测试方法为测定某些特 殊参数提供了最准确或最实用的途径；而在某些情况下，甚至成为唯一可靠的 调查技术。 借鉴现场调查和物理实验中示踪剂的概念，本文假定在泥沙颗粒上吸附某 种不降解的，完全跟随泥沙运动的物质，以此作为一种理想的示踪剂，来研究 第一章绪论 沉积物的迁移扩散情况。 1 3本文主要工作 本文属于国家自然科学基金的一部分，最终的目标是研究吸附于沉积物中的 污染物的迁移转化问题。而本文主要开展了以下几个方面的工作： 1 建立了s w a n 风浪模型及二维潮流模型，并用经典的实验结果和渤海湾 的实测数据对模型进行了验证。 2 建立了二维对流扩散及悬沙模型，分别模拟了在潮流单独作用下以及波 浪和潮流共同作用下的渤海湾泥沙运动情况，并用实测资料进行了验证。 3 模拟了渤海湾的流场和余流场，初步分析了渤海湾沉积物可能的迁移扩 散情况；进一步模拟了保守物质在渤海湾的迁移扩散情况，并与泥沙的运动进行 了定性的比较分析。 4 建立了沉积物示踪模型，用该模型模拟了天津近岸海域沉积物的迁移扩 散情况，并与根据余流场分析的渤海湾沉积物可能的迁移扩散情况、保守物质迁 移扩散情况及泥沙的运动情况进行了对比分析，初步得出了沉积物迁移扩散的一 些特性。 第二章数学模型的原理及数值计算方法 第二章数学模型的原理及数值计算方法 河口海岸沉积物运动是风、浪、流共同作用的结果，也是海水运动和海岸工 程建筑物及海岸带相互作用结果。研究水体中污染物的运动和沉积物的运动需要 建立合理的风浪模型、潮流模型、对流扩散模型和泥沙模型。而要跟踪沉积物的 迁移扩散，则需采用一定的示踪技术。本章建立了风浪模型、潮流模型、对流扩 散模型、泥沙模型以及沉积物示踪模型。 2 1风浪模型 2 1 1 s w a n 数学模型简介 在基于动谱平衡方程的第三代波浪数值模型中，s w a n 模型是应用较广、较为 成熟的模型之一。s w a n 模型适用于深水、过渡水深和浅水多种水深情况；模型 包括能量输入、损耗和非线性相互作用机理，源项的处理应用当今海浪研究最新 成果，尤其在非线性项中加入三相波相互作用项，能合理模拟近岸波浪传播的周 期变化；将随机波浪以不规定谱型的方向谱表示，更接近真实海浪；模型计算的 边界条件比较简单，不要求闭合边界条件，只要适当选择计算域的边界，即便只 能确定其中一个迎浪边界条件，也能获得可靠的效果。 但是s w a n 模型也有其局限性，即：模型不包含绕射，不适于模拟绕射效应显 著的港湾波浪。对于不规则短峰波，现场观测及试验均证明，除障碍物顶端后一 至两个波长范围内本模型模拟结果有所变小外，当入射方向谱不太窄时，绕射效 应是比较小的，因此模拟大面积的风浪情况时，此缺点并不十分明显。 在对渤海湾这种以风浪为主的大范围浅水波浪场进行模拟时，s w a n 模型是合 理的选择，能够用于各种复杂的气象和地形条件。 2 1 2 s w a n 模型的数学表示 2 1 2 1波浪的传播方程 作用量密度的变化率可以用作用量平衡方程吲表示： 妄+ 丢e + 专c + 昙c ，+ 品g = 軎 c 2 第二章数学模型的原理及数值计算方法 其中：方程左边第一项表示作用量密度随时间的变化率。第二项和第三项表示作 用量密度在几何空间的传播( 变化率c 。和c ，) 。第四项表示由流和变化的水深弓i 起的折射和变浅作用( 传播速度q ) 。第五项表示流和变化的水深引起的频移( 传 播速度c 占) 。其中n ( t r ，口) 为动谱密度( 能谱密度e ( t r ，臼) 与相对频率盯之比拍们) ； 仃为波浪的相对频率( 随水流运动的坐标系中观测到的频率) ；秒为波向( 各谱 分量中垂直于波峰线的方向) ；c x ，c y ，巳，c e 分别表示x ，y 方向以及仃， 秒空间的波浪传播速度。方程右边的s 表示能量源项，它可写成几个不同类型的 源项之和s ，= s 加+ 黾+ s 其中s 折表示风输入项；s 出表示由白浪破碎、底摩 擦、变浅破碎引起的耗散作用项；s 。，表示波和波之间的非线性相互作用项。 从线性波动的波包理论和波峰守恒率可以在几何空间和谱空间得出下面的 变化率m 1 ： 鲁= ( ) = 西肌扑商赫睁痧 ( 2 - 2 ) 百d o = g = 塑a di l 型o t + 疗v ，dl _ c g k 詈 ( 2 - 3 ) 一d o g 一一- 1l o 甜t r 丽a d 蒋t - 筹i 刚2-4dtk ) 一= c 。= 一- 一庀ii) 9 a d 锄a 惭l 、7 其中：d 是水深，m 是仃方向的坐标轴，s 是秒方向的坐标轴。另外随体导数表 示为瓦d = 昙+ ( 西+ 西) v - 其中v - = ( 昙，杀) ，未和昙分别代表求仃和口方向 的方向导数，云代表波数单位矢量。 2 1 2 2 源项与汇项的输入 ( 1 )风能输入 风能输入项根据c a v a l e r i & m a l a n o t t e r i z z o l i 、k o m e n e t a l m l 的研究成果得出： 风输入项是线性增长和指数增长之和： s i n ( t r ，口) = a + b e ( o ，秒)( 2 - 5 ) 其中：a 表示线性增长，b e 表示指数增长。 ( 2 )波能耗散 波能的耗散机制分为三种：白浪破碎、底摩擦和变浅破碎。在深水情况，风 浪的白浪破碎起主要作用，控制着谱的高频部分的饱和程度；在中等深度和浅水 情况下，底摩擦变得重要，成为波能耗散的主要因素；在浅水破碎带附近时，因 变浅引起的波浪破碎起主要作用。 ( 3 ) 非线性波波相互作用 非线性相互作用后得结果是使能量在谱内重新分配。在深水情况四相波相互 第二章数学模型的原理及数值计算方法 作用比较重要，浅水情况三相波相互作用变得重要。其中三相波相互作用引起的 能量交换采用r i s e t a l 的研究成果；四相波相互作用的处理采用h a s s e l m a n n 的离 散迭代近似法( d i a d i s c r e t ei n t e r a c t i o na p p r o x i m a t i o n ) 。 2 1 &计算域边界条件的处理 计算域边界可以分为陆地边界和水域边界。陆地边界认为能将入射波吸收而 不产生波浪反射；迎浪面边界条件一般根据现场观测或波浪模型数值模拟得到。 通常现场观测到波浪数据较少而且不一定在迎浪边界上，此时通过其它大尺度波 浪数学模型能得到可以接受的误差范围内边界波浪数据。 2 2潮流模型 自然界中的水体运动，实质上都是三维问题，但是为了简化数值模拟和节约 计算工作量，可以根据其物理特性的不同分为一维、二维或三维问题来研究。本 文的研究对象为水深较浅的河口、近岸、海湾等处的水体，其垂向尺度远小于水 平尺度，流速在垂直方向的大小和变化都远小于水平方向的大小和变化，其流动 特征可用沿垂向积分后的物理量表示。对于水深较小近岸、海湾处的水体其水平 尺度( x ，y 方向) 一般在1 0 3 1 0 5 m 量阶，而其垂向( z 方向) 尺度一般在 1 0 1 1 0 2 m 量阶，可将其水体运动视为浅水流动。故本文针对渤海和渤海湾，应 用二维水深积分的长波方程来研究海水的运动。 2 2 1潮流模型的控制方程 沿水深积分的二维潮流模型的控制方程为如f 偏微分方程组： 连续性方程质量守恒方程： 譬+ 罢+ 粤喊 x 方向动量( 守恒) 方程： 至+ a ( , a p u ) + a ( 蜀p v ) ：f q 一旦堡一堕+ 垒一鱼 a t i g x a ypi g x 。 i g x pp + 吖窘寺卜 y 方向动量( 守恒) 方程： 8 ( 2 - 6 a ) f 2 6 b ) 第二章数学模型的原理及数值计算方法 一a q + 旦型! + o ( p q v ) ：一廖一一h 堡一g h 丝- i 鱼一旦 研， 出 、o y po y砂pp ( 2 6 c ) 一一 、 、一， 托( 罄+ 罄卜最 其中：x 、y 为水平方向坐标，f 表示时问，f ( x ，y ，t ) 为水位，u 、矿分别为沿 水深平均的流速在x 、y 方向的分量，p 、q 分别为流体在“y 方向的单宽流 量，h ( x ，y ，t ) = f ( x ，y ，t ) + h ( x ，y ，t ) 为水深，h ( x ，y ，t ) 为相对于基准面的水深，s 。、 u 。、圪分别为单位水平面积上的源项强度及该源项初始速度在x 、y 方向的分 量，卢为非均匀流的动量修正系数( 由流速垂向分布引起) ，厂= 2 0 s i n 0 为柯氏 力( c o r i o l i s ) 系数( 由地球自传引起) ，国为地球自转角频率，0 为纬度，g 为 重力加速度，本文取g = 9 8 1 m s 2 ，夕为水体密度，本文取p = 1 0 2 3 k g m 3 ，尼为 自由水面处大气压强，f 。、f 。分别表示自由水面剪切应力( 风应力) 在x 、y 方 向的分量，、分别表示水体底部剪切应力( 底摩擦应力) 在工、y 方向的 分量，占为水深平均涡粘系数。 2 2 2 潮流模型初始条件与边界条件 对于水动力学模型，在模拟初始时刻( f = 0 ) 应给定水位、单宽流量的初始值， 即： f ( x ，y ，0 ) = f o ( x ，y ) ( 2 7 a ) p ( x ，y ，0 ) = p o ( 工，y ) ( 2 - 7 b ) q ( x ，y ，0 ) = q o ( x ，少) ( 2 7 c ) 其中：上标( o ) 表示初始时刻已知的全场数值。模型的启动方式可根据给定初 始时刻水位或( 单宽) 流量，启动方式分为：冷启动和热启动。冷启动是指：初 始时刻全场的单宽流量取零，水位取为统一值( 最高水位、最低水位或平均水位) 。 热启动是指：初始时刻全场的水位、( 单宽) 流量由试算或其它计算结果给出。 如无特殊说明本文均采用冷启动的方法，水位取平均水位。 水动力学模型在模拟计算过程中需要指定其边界条件。边界条件分为：闭边 界( r o ) 条件和开边界( r 1 ) 条件。在闭边界处流速采用无渗透、部分滑移边 界条件。无渗透是指：沿边界法线方向( 元l ) 流速为零，即： 圪| r n = v n v 0 = 0 ( 2 8 ) 部分滑移是指：闭边界处沿边界切线方向( 气，) 存在切向流速即旷气0 。 其切向流速沿法向的导数为： 刮r = v 妒气) 民诋 ( 2 - 9 a ) 第二章数学模型的原理及数值计算方法 其中：口r 0 【0 ，1 为滑移系数，其取值情况反映了边界处流速的滑移状态，有如 下关系： j 口r 0 = 1 滑移边界条件 0 口r 0 o 5 m m ( 2 2 1 ) 【v 5 0驯 5 0 匕，一= 两t w ( 2 2 2 ) p 为( 蔫+ 成篙 弘2 4 ， 第二章数学模型的原理及数值计算方法 其中：刀为曼宁系数( 适用范围0 0 1 2 - - 0 0 2 5 ) 礅口o = 0 0 2 3 ；徽, o o = o 0 4 l ， 其中l 为波浪底摩阻系数，适用取值范围0 0 0 7 5 - 0 0 1 5 ，近似取厂w = 0 0 1 可得 成= 0 0 0 0 4 。f f 式( 2 2 4 ) 适用的含沙量范围：0 0 3 k g m 3 4 0 k g m 3 。 2 4 3河床变形方程 由于水流动力条件变化或泥沙边界条件的变化会引起河床的冲刷和淤积。随 着涨落潮流速的变化，河床一直处于不断调整以适应水动力条件和边界条件的动 态变化之中。因此，研究河床的变化情况也具有非常重要的意义。河床变形方程 同样可以用质量守恒方程来推导。方程描述的物理意义为：单位时间进入微小水 柱内含沙量的变化和河床上泥沙的淤积量之和。方程形式为： 孕+ 士f 监+ 盟1 ：0 ( 2 - 2 5 ) 祝 l 一lo xa y 其中：p 。、q ，( m 2 s ) 分别代表石和y 方向的总体积输沙率；刀，( 无量纲) 为 底床孔隙率；乙为底部高程。 在忽略推移质运动的假定下： p 。：型( 2 2 o ) s = 三一 【z ps g ，：笙 ( 2 2 7 ) g ，= 二一 u 。，) ps 其中：肛( k g m 3 ) 为泥沙密度；p 、g ( m 2 s ) 为x 和y 方向的单宽流量；s ( k g m 3 ) 为水深平均的悬沙质量浓度。 2 5沉积物示踪模型 在现实世界中，吸附于泥沙颗粒上的放射性同位素、难降解的污染物( 如重 金属) 等均可在一定程度上看作沉积物运动的示踪剂。本文借鉴现实示踪剂的概 念，提出了用于数值模拟的理想的示踪剂，该理想示踪剂满足以下条件： 1 ) 完全吸附于泥沙颗粒上，即只存在于固相中，而不存在于水相中； 2 ) 没有质量和大小，完全跟随泥沙颗粒运动； 3 ) 为保守物质，即不降解。 第二章数学模型的原理及数值计算方法 掣+ 掣+ 掣i b 鲁一x y 引 陆2 8 ， 国叙 咖缸i “反咖j 一杀旧i a c p + h d 叫等卜 一 醚：q d e v pdeposii。n(2-29)。 【q e r o g e r o s i o n 2 6数值求解方法 2 6 1s w a n 模型的离散格式与求解 s w a n 模型采用了全隐式有限差分格式，无条件稳定，即使在很浅水域，也 可以取较大时间步长。时间以简单的常数时间步长垃离散化。这与w a m 模型 和w a v e w a t c hi i i 模型的时间离散化是不同的，在这两个模型中，波传播的时 间步长与源项的时间步长是不同的。地理空间使用矩形网格，在x 方向和y 方向 风别取常量步长缸和衄。模型中的谱以常量a 8 离散化( 对数率分布) 。可以选择 只在先定义好的方向因数内计算波结构( 氏；。 秒 氏。，即在限定方向内面向海 岸传播) 。离散频率定义在一个确定的低频分离点和一个确定的高频分离点之间 ( 谱的预知部分) 。在低频分离点以下( 如开阔区域条件下名；。= 0 0 4 h z ) 假设谱密 度为0 。在高频分离点以上( 如开阔区域条件下1h z ) 加入诊断尾迹厂”( 用于计 算高频非线性波一波相互作用以及计算积分波参数) 。应该在4 和5 之间取值。如 果使用k o m e n 的风输入方程则m = 4 ，j a n s s e n 的风输入方程则m = 5 。 动谱平衡方程进行离散后的差分方程表示为： 第二章数学模型的原理及数值计算方法 j 。j ，j e + i n l x j 4 如 + 阻坚垃娑堕竺生划“”+ ( 2 - 3 2 ) l2 仃i 、7 f二型坚兰垒生生l二兰丛2垒a暑掣立掣：二矗=l剖2：i：，岛 口ll 仃l ；，； 上式中：f ，为时间层标号5 t ，f ，如为工，y ，盯和秒方向响应网格标号5 出为时间步长； 缸，衄风别为为地理空间上x 、y 方向的空间步长； a t r ，a 8 为谱空间相对频率仃和方向分布p 的步长5 n 为每时间层迭代次数。 方程右边源汇项中的行为刀或玎一1 ；r 【0 ，1 】，y 【0 ，1 】，系数y 和7 7 的取值 大小决定了谱空间的差分格式是偏于迎风格式还是偏于中心格式，因此决定了在 谱空间和方向空间的数值离散程度和收敛性强弱。如当y = 0 或卵= 0 时为中心差 分格式，数值离散趋于0 ，计算准确度最高；当l ，= 1 或玎= 1 时为迎风差分格式， 数值离散程度最大，但收敛性最好。 2 6 2 潮流模型的离散格式与求解 本文采用空间交错网格系统对空间进行离散，将空间划分成许多单元，其中 水位节点和流速( 或单宽流量) 节点交错布置于单元的中心和两边，地形节点位 于单元的角点，其示意图如图2 1 所示。采用空间交错网格系统可以有效的避免 求解过程中数值振荡的出现( f l e t c h e r , c a ，1 9 9 1 ) 。将地形节点布置于网格角点 可以在固定网格空间步长情况下描述最详细的地形情况，其他位置的地形由角点 处地形插值而得。 用a d i 法( a l t e r n a t i n gd i r e c t i o ni m p l i c i tt e c h n i q u e ，即交替方向隐格式法) 对潮流模型及对流扩散模型的控制方程进行求解。a d i 法的特点是将每一时间步 长分为两个半步长，将一个二维问题转换成多次求解一维问题，从而避免了求解 二维问题的全元素矩阵。结合边界条件，半时间步长内差分方程可以化为三对角 方程组，可采用追赶法求解，避免了直接求解大型稀疏矩阵。a d i 法具有无条件 稳定、计算速度快、存储量小和适用范围广的优点。 第二章数学模型的原理及数值计算方法 吉凄 s2y 口+ 白 嘛瑚：、q ；绁lc 乙 、叠 12 l y 团t b - 噙l ，2 二l b 刷j 黜l j 乙 吉l 帆艘s2 丫 口岁曲 基叁 ，卜l ，2 “3 ，工j 图2 - 1 空间交错网格系统 ：地形节点；口：水位点： o ：流量、流速节点；_ ：流量方向； 0 ，1 ，2 - ( 地形) 节点标号，单元编号 前半个时间步长( n 专n + 1 2 ) 在前半个时间步长内，将连续性方程与工方向动量方程联立，对f ，p 进行 隐式求解，对g 进行显式求解。 连续性方程离散为： & | 2斌学镀 ( 2 3 3 ) 其中：血，缈，垃分别为x ，y 方向空间步长和时间步长，下标i ，j 分别代表x ，y 方 向单元编号，上标n 代表时间步长数，即，z ，n + 1 2 分别代表时刻 t = n a t ，t = ( n + 1 2 ) a t 。整理上式可得： 口2 ，- l p f n 一+ l l 2 2 ，+ 6 2 f l f 孑“2 + c 2 f l p f n + l + l 2 2 ，= d 2 f l ( 2 3 4 ) 其中： 口z 一石 b 2 i l = 2 a t c ：h2 石 - 2 劬一崇q 7 , ：+ 1 2 - q ) + 出乜 ( 2 - 3 5 a ) ( 2 3 5 b ) ( 2 - 3 5 c ) ( 2 3 5 d ) 第二章数学模型的原理及数值计算方法 对于动量方程中的非线性项，即使在应用隐格式的情况下，也会产生不稳定 现象( w e r e ，1 9 7 6 ) 。消除其不稳定性的方法有，应用三个时间层对进行时间中心 差( f a l c o n e r , 1 9 7 7 ) 或对流速进行光滑处理( p o n c e ，1 9 8 1 ) 。另外一个方法是对非线 性项和涡粘项进行时间中心迭代处理( r o a c h e ，1 9 7 6 ) ，而对流交叉项应用一阶迎 风格式，这样可以消除流速梯度较大区域内的局部振荡现象( o w e n s ，1 9 8 7 ) 。本文 x 方向动量方程离散为( f a l c o n e r 等，1 9 9 9 ) ： p ：墓| ；j ? 一p n i “- l | 7 l 2 。i = 厂，鼓l ，2 ， h 三t ，2 ，( 名) 知，一( 尼) ! ， p a x ( p 矿) 0 l 2 ，+ l 2 一( p y ) ：l ，2 ，一l 2 少 一g - - , + l ，2 ，j( 盟笋 一( 罕 仁3 6 ， 一 g丝n+l学2 n - l7 2 雁五万丽 + f “墨l ，2 1 j l 辰z 五i i ：z 了。r 一 + 占( 亟鼍竽亟+ 亟生产 兵中：上标( ) 代表迭代修正僵例如，上标( i ) 代表周围网格点的半均僵 上标( ；) 代表依照迎风格式取值，举例表示如下： f p ”“2 第一次迭代 2 1 丢g 肘z + p n - l 2 ) 以后各迭代( 2 - 3 7 ) 霸位川，：= 昙眈知：+ 川，：) ( 2 3 8 ) 扎- - n m 广爿华+ 学】 亿3 9 ， p i + 1 1 2 , i = 脍鼍魏猫0 弘4 。， 1 p 川当2 ， 峥7 整理动量差分方程可得： 口2 f f 孑2 + 6 2 ，p 州n + l ，2 ，2 + c 2 ，嚣y 2 = d 2 f ( 2 - 4 1 ) 其中： 口2 f = 一g h _ a t ( 2 4 2 a ) 第二章数学模型的原理及数值计算方法 咖i n 2-n 2一 = 尝 九= i n2 ( - - n 、2 ：k ， 一俘-i2 a t l = - ：_ 熹w 靴； 娥m ，( 名) 一( 名) ：， p 厶x + ( _ 。) ：l ，2 ， 一石, tb l n 垆l 一何m ) + 2 氍鬟麓h ； ( 2 - 4 2 b ) ( 2 4 2 c ) ( 2 4 2 d ) 口2 ( f 1 ) n + 1 1 2 + 如( f _ 1 ) p f - 1 n + l ，2 1 2 + c 2 ( ) 甜，2 = d 2 ( “) 口2 f l p f n - l + i 2 2 ，+ 也h筝才，2 +