预估修正法对河口海岸海洋模式稳定性的提高3

预估修正法对河口海岸海洋模式稳定性的提高3朱建荣 杨陇慧 朱首贤华东师范大学河口海岸国家重点实验室 上海提要 河口海岸海洋模式时间差分为欧拉前差格式 在涡动粘滞系数较小的情况下 模式存在着弱不稳定性∀为提高模式的稳定性 采用预估修正法对模式中科氏力项作半隐式处理∀首先在理论上证明了模式稳定性的提高 随后设计一个高分辨率网格 考虑实际岸线和水深!边界通量!密度梯度力和风应力 把模式应用于冬季渤!黄!东海环流的数值模拟 数值试验也证明了模式稳定性的提高∀改进后的河口海岸海洋模式较成功地模拟出了冬季渤!黄!东海环流 黑潮!台湾暖流!对马暖流!黄海暖流和沿岸流 与实际观测和研究结果较为一致∀关键词 河口海岸海洋模式 预估修正法 模式稳定性 渤!黄!东海环流中图分类号 ° 1河口海岸海洋模式 ∞≤ 2 ∞ ∏ ≤ ≥ 2≥ 是当今国内外应用较为广泛的海洋模式 ∏ εταλ 许朋柱等 但在涡动粘滞系数较小的情况 尤其在无潮汐混合和层结存在的情况下 ∞≤ 2存在着弱不稳定性∀理论上可证明在无粘的情况下 欧拉前差格式是惯性不稳定的∀本文在原∞≤ 2 的基础上 采用预估修正法 ° 2≤ ≥ 对模式中科氏力项作半隐式处理 以提高模式稳定性∀把改进后的模式应用于冬季渤!黄!东海环流的数值模拟 理论分析和数值试验均表明了预估修正法提高了∞≤ 2的稳定性∀以往对长江河口作了许多工作 赵世清 于克俊 许朋柱等 刘兴泉 改进后的模式达到了提高稳定性的预期目的 为长江河口动力过程研究提供了先进的研究手段∀1 ΕΧΟΜ−σι和预估修正法111 非正交坐标系下河口海岸海洋控制方程组引入水平非正交曲线 Ν Γ 和垂向Ρ坐标系 曲线坐标 Ν Γ Ρ 定义为Ν Ν ξ ψ Γ Γ ξ ψ Ρζ ΦΗ Φ其中 ξ!ψ和ζ为向东!向北和向上的笛卡儿坐标轴 垂向坐标Ρ从海底 ζ η 变至海表面 ζ Φ Φ为海表面波动 Η为总水深∀河口海岸海洋控制方程组包括动量!质量连续!温度!盐度和密度过程∀Ν Γ方向动量方程为3国家自然科学基金资助项目 号 高等学校博士学科点专项科研基金资助项目 号∀朱建荣男 出生于 年 月 博士 教授 ∞2 ∏ 1 ∏1 ∏1收稿日期 2 2 收修改稿日期 2 2第 卷第 期年 月海 洋 与 湖 沼≤∞∞×≥ ≤∂5∆ϑυ 5τ5∆ϑΥ υ5Ν5∆ϑς υ5Γ5ϑΞυ5Ρ∆η ς ϖ55Νϑηυ55Γϑηϑφ ∆η υ ϖ55Νηη ηη γ∆5Φ5Νγη ∆Θο5∆5ΝΘΡΡ5Θ5ΡδΡγη ∆Θο55ΝΘΡΘδΡ∆55ΡΚµ5ϑυ5Ρ∆ϑΦυ5∆ϑϖ 5τ5∆ϑΥ ϖ5Ν5∆ϑς ϖ5Γ5ϑΞϖ5Ρ∆η Υ ϖ55Νϑηυ55Γϑηϑφ ∆η υ ϖ55Γηη ηη γ∆5Φ5Γγη ∆Θο5∆5ΓΘΡΡ5Θ5ΡδΡγη ∆Θο55ΓΘΡΘδΡ∆55ΡΚµ5ϑϖ5Ρ∆ϑΦϖ方程 ! 中 Ν和Γ方向速度分量 定义为υϖ 可表示为 υηϑξΝυψΝϖ ϖ ηϑξΓυ ψΓϖ 其中 υ!ϖ为ξ!ψ方向速度分量 ΝξψΓϑΝψξΓϑΓξψΝϑΓψξΝϑϑ为雅可比函数 表示为 ϑ ξΝψΓ ξΓψΝ 下标Ν!Γ表示微分运算 坐标变换后的因子η!η 定义为η ξ Ν ψ Ν η ξ Γ ψ Γ Υϑη υηηϖ ςϑ η ϖηη υ其中 η ψΝψΓ ξΝξΓ φ为科氏参数 γ为重力加速度 Κµ为垂向涡动粘滞系数 Φυ!Φϖ分别代表水平动量扩散项 Θ和Θο为扰动和参考密度 它们满足Θ Θ Θο关系∀质量连续方程为 5Φ5τϑ55Ν∆ϑΥ55Γ∆ϑς5Ξ5Ρ温度!盐度方程为 5ϑ∆Η5τ5ϑ∆Υ Η5Ν5ϑ∆ς Η5Γ5ϑΞΗ5Ρ∆55ΡΚη5ϑΗ5Ρ∆ϑΦΗ5ϑ∆σ5τ5ϑ∆Υ σ5Ν5ϑ∆ς σ5Γ5ϑΞσ5Ρ∆55ΡΚη5ϑσ5Ρ∆ϑΦσ式中 Η为位温 σ为盐度 Κη为热力垂向涡动扩散系数 ΦΗ!Φσ分别代表热量和盐度水平扩散项∀密度方程为 ΘΘ Η σ其中 Ρ坐标中垂向速度Ξ为 Ξ ω ΡΥ 5∆5Νς5∆5ΓΡ5Φ5τΥ5Φ5Νς5Φ5ΓΦυ!Φϖ!ΦΗ和Φσ中水平扩散系数由≥ 公式 计算 水平网格元大小和流的水平切变成正比 Κµ和Κη使用改进后 等 ! 1 阶湍流闭合模型计算∀在本文中 模式仅考虑热力过程对流场的影响 而不考虑流场对温盐的影响∀海表面和海底动量和热力边界条件的给出可参见 ∏ 等 ∀112预估修正法在∞≤ 2 中 动量方程分以下三步做时间积分∀第一步∆ϑυνα ∆ϑυν∃τ5∆ϑΥ νυν5Ν5∆ϑς νυν5Γ5ϑΞνυν5Ρ∆η ς νϖν55Νϑη海洋与湖沼 卷υν55Γϑη ∆η ςνϑφ ∆ηυν ϖν55Νηη ηγη ∆Θο5∆5ΝΘΡΡ5Θν5ΡδΡγη ∆Θο55ΝΘΡΘνδΡ ∆ϑΦνυ∆ϑϖνα ∆ϑϖν∃τ5∆ϑΥ νϖν5Ν5∆ϑς νϖν5Γ5ϑΞνϖν5Ρ∆η Υ νϖν55Νϑηυν55Γϑη ∆η Υνϑφ ∆ηυν ϖν55Γηη ηγη ∆Θο5∆5ΓΘΡΡ5Θν5ΡδΡγη ∆Θο55ΓΘΡΘνδΡ ∆ϑΦνϖ 第二步∆ϑυνβ ∆ϑυνα∃τ∆55ΡΚµ5ϑυνβ5Ρ∆ϑϖνβ ∆ϑϖνα∃τ∆55ΡΚµ5ϑϖνβ5Ρ第三步∆ϑυν∆ϑυνβ∃τη γ∆5Φν5Ν∆ϑϖν∆ϑϖνβ∃τη γ∆5Φν5Γ在第一步中 仅考虑非线性项!科氏力项!密度梯度力项!水平涡动粘滞项和由坐标变换产生的曲率项 所有这些项均作显式处理∀可以证明 见附录Α 欧拉前差和科氏力为显式格式 在无粘或粘滞系数较小的情况下 存在着惯性弱不稳定性 • εταλ ∀第二步中 仅考虑垂向涡动粘滞项 并做隐式处理 这样在垂直方向可及大地提高分辨率∀第三步中 仅考虑外重力波产生的水位梯度力 这个快过程显然作隐式处理 把υν!ϖν代入连续方程 用半隐格式求解水位分布 消除了≤ƒ 判据的严格限制 提高了计算效率∀求出水位分布后 由第三步可直接得出三维速度场的分布∀在不改变∞≤ 2 中二步时间层算法的基础上 引入预估修正法改进第一步的计算格式 可提高模式的稳定性 见附录Β证明 ∀预估时间步为∆ϑυ3να ∆ϑυν∃τ∆η ς νϑφ ΡΗΤν∆ϑϖ3να ∆ϑϖν∃τ∆η Υ νϑφ ΡΗΤν其中 ΡΗΤ为第一步中右边除科氏力项外所有其他项∀求出υ3να ϖ3να后 继续上述第二!三步 可得到预估步的速度场υ3νϖ3ν∀修正步为∆ϑυνα ∆ϑυν∃τ ∆η ≈Βς 3ν Β ς ν ϑφ ΡΗΤν期朱建荣等 预估修正法对河口海岸海洋模式稳定性的提高∆ϑϖνα ∆ϑϖν∃τ∆η ≈ΒΥ 3ν Β Υ ν ϑφ ΡΗΤν取Β 1 φ∃τ可保持惯性波的中性放大∀求出修正步中第一步的速度值υν αϖν α后 继续上述第二!三步 最后可获得一个时间层次上的水位Φν 和速度υνϖν分布∀2数值试验211未作预估修正取Β 即为未作预估修正的原非正交曲线坐标系下∞≤ 2 应用于冬季渤!黄!东海环流的数值试验∀模式区域包括整个渤!黄!东海 考虑实际海岸线和水深 设计高分辨率的水平网格 长江口外分辨率为 ) 其余区域为 ) 图略 ∀垂向分 层 时间步长取 ∀初始流场和水位取为υϖ 和Φ ∀侧向固边界取法向速度为零 台湾海峡!台湾以东海区!吐噶喇海峡!对马海峡和长江作水侧边界处理 其余海区作固边界处理∀基于以往的观测和研究结果 ƒ εταλ ≠∏ εταλ εταλ εταλ 台湾海峡!台湾以东海区!吐噶喇海峡!对马海峡和长江的开边界条件的通量分别取为 1 ! 1 ! 1 和 1 ! 1 Σϖ Σϖ 流进和流出的流量总体保持平衡∀驱动环流的动力因子为边界通量!密度梯度力和风应力∀模式的初始温度和盐度场取自渤海!黄海!东海海洋水文图集 海洋图集编委会 月各标准层上 ! ! ! ! ! ! 及底层 温度和盐度分布 再水平插值到模式网格点和垂向插值到各Ρ面上∀本文不考虑温度和盐度随时间的变化 只考虑它们对流场的影响∀风应力参数取自渤海!黄海!东海海洋气象图集 海洋图集编委会 月资料 再转化为风应力值 εταλ ∀时间步长取 ∀运行未作预估修正的∞≤ 2 计算到 模式天时 黑潮主轴及附近开始出现乱流现象 不稳定性 图略 ∀212已作预估修正采用上述未作预估修正法时数值试验相同的初始和边界条件 将已作预估修正后的∞≤ 2 应用于冬季渤!黄!东海环流的数值试验∀采用预估修正后 模式的时间步长还是受一定限制的∀由附录Β中的证明可知 稳定性条件必须同时满足Β 1 φ∃τφ∃τ ∀由φ∃τ ∃τ φΥ ∗ 这个条件一般是容易满足的∀取Β 1 为和未作预估修正法时的试验结果比较 时间步长仍取 ∀运行作过预估修正的∞≤ 2 计算到 模式天时 没有出现不稳定现象 继续运行模式 模式天 模式仍然保持稳定∀213比较以上两个数值试验表明预估修正法提高了∞≤ 2 稳定性∀图 ! 分别为模式运行至 模式天时 渤!黄!东海表!底层流场分布∀模式较好地再现了黑潮!台湾暖流!对马暖流!黄海暖流和各沿岸流∀近岸浅海表层流主要被风生流所控制 闽浙沿岸流在东北风的作用下 沿岸向南流动∀向南流动的西朝鲜沿岸流在离开西朝鲜海岸后 流入 海洋与湖沼 卷九州海峡 随后并入对马暖流∀表层黑潮沿陆架坡朝东北方向流动 流速大约为 Ù ∀在底层 近岸浅海海流较弱 在陆架坡处 黑潮垂向流速急速减少 在陆架坡外侧甚至流向相反 这表明黑潮具有很强的斜压性∀从表!底层均可看到 对马暖流来自台湾暖流 而不是来自九州西南黑潮的分离 因为在琉球深沟的北部有一小流速区∀这个结论和 等 描述的东海和黄海的环流形态是一致的∀台湾暖流来自台湾海峡和台湾东北侧 且主要来自台湾东北黑潮的入侵∀ 等 曾根据大量卫星跟踪浮标轨迹得出了冬季表层东海!黄海东部流场分布 将本文模拟结果与之比较 可见流系的位置和量值与 等 卫星调查结果均较为相似∀图 模式模拟的 月渤!黄!东海表层 和底层 流场ƒ 1 ≥ ∏ ∏ ∏!≠ ∞ ≤ ≥ ƒ ∏3 结语本文研究在理论上证明了∞≤ 2在采用时间差分为欧拉前差!科氏力为显式格式!涡动粘滞系数较小的情况下 存在着弱不稳定性∀采用预估修正法对模式中科氏力项作半隐式处理 先在理论上证明了模式稳定性的提高∀又设计一个高分辨率网格 考虑实际岸线和水深 综合考虑边界通量!密度梯度力和风应力 把模式应用于冬季渤!黄!东海环流的数值试验∀实验结果表明 预估修正法提高了模式的稳定性∀改进后的∞2≤ 2较成功地模拟出了冬季渤!黄!东海环流 黑潮!台湾暖流!对马暖流!黄海暖流和沿岸流 与以往实际观测和研究结果较为一致∀期朱建荣等 预估修正法对河口海岸海洋模式稳定性的提高附录Α欧拉前差惯性弱不稳定性分析为简化问题 突出科氏力项的不同有限差分格式 只考虑线性动量方程 雷诺应力项也考虑简单的形式5υ5τφϖ ρυ5ϖ5τφυ ρϖ 写成复数的形式则为5ω5τιφ ρ ω 其中 ω υ ιϖ ι 其余变量如常用∀对 作二步时间层的欧拉前差 并仅对科氏力作半隐处理 可得ων Ρ ων ιΦ≈Βων Β ων 其中 [Β[ Φ φ∃τ Ρ ρ∃τ ν和ν 分别表前一和当前时间步∀若Β 则为欧拉前差 显式格式 若Β 1 则为半隐格式 若Β 则为欧拉后差 隐式格式 ∀为了得到数值稳定性判据 假设ων Κων 可见 Κ 数值解是中性的 Κ 数值解是衰减的 Κ 数值解是增长 即不稳定的∀在∞≤ 中 科氏力是作显式处理的 Β 这里仅讨论Β 情形∀把 代入Κ Ρ ιΦ取Κ绝对值 ÞΚÞΦΡ Φ 如ρ如ρΞ由 可见在无摩擦 ρ 情况下 数值解是无条件不稳定的∀若ρΞ 且Φ远小于 Φφ∃τ ∃τφΥ ∗ 这个条件一般是容易满足的 可推算出稳定性条件Φ Ù [Ρ[ Φ Ù即数值模式的耗散必须足够大才能抑止由于科氏力项不恰当的差分引起的数值不稳定∀在无潮汐混合和层结存在抑止垂向混合情况下 ∞≤ 中垂向和侧向湍流涡动系数可变得很小 会出现弱不稳定的情形∀附录Β预估修正法对模式稳定性的提高对方程 分别作预估步和修正步差分如下ω3ν Ρ ων ιΦωνων Ρ ων ιΦ≈Βω3ν Β ων 右上角带3值表示预估值 [Β[ 将预估步代入修正步可得ων ≈ Ρ ΒΦ ιΦ ΡΒ ων 将 代入 可得到特征值Κ关系式 Κ Ρ ΒΦ ιΦ ΒΡ 对 取绝对值 可得 ÞΚÞ Ρ Φ Β Φ ΒΡ 现考虑无摩擦极端情况 ρÞΚÞ ΦΦΒ Φ Β ΦΦ Φ若Β若Β若ΒΕ以及Φ[若Β 以及Φ[海洋与湖沼 卷其中 Ε为一小正数∀若Β 格式是无条件不稳定的 和欧拉前差相似 若Β 1 格式仍是不稳定的 为保证格式稳定 令ΒΕ 代入 且设ÞΚÞ 在假设Φ[ 的条件下 如前所述 这个条件是容易满足的 可得到小正数Ε Φ Ù Φ∀也就是说 采用预估修正法后 即使流体不存在摩擦 只要取Β 1 Φ ÙΦ这个格式是稳定的∀这从理论上证明了与欧拉前差格式相比较 预估修正法提高了模式的稳定性∀若Β 以及Φ[ 模式也是稳定的∀若ρΞ 即存在摩擦的情况下 与ρ 情况比较 预估修正法将更趋稳定∀参 考 文 献于克俊 1长江口余流和盐度的三维数值计算1海洋与湖沼 ) 许朋柱 毛锐 1长江口南支潮流的一个数值模拟1海洋与湖沼 )刘兴泉 1沿岸海区冬季垂直环流及其温盐结构的数值研究 ¬1环流基本特征1海洋与湖沼 ) 赵世清 1长江口潮流的一种数值模式1海洋与湖沼 )海洋图集编委会 1渤海!黄海!东海海洋图集 水文1北京 海洋出版社 ) 海洋图集编委会 渤海!黄海!东海海洋图集 气象1北京 海洋出版社 )≤ ∏ ≠∏ εταλ 1⁄ ≤ ≠ √ ∞ ≤≥ 1≤ ≥ Ù) ∏ ƒ 1 2 ∏× 2 • ⁄≤ )ƒ ∏≠ 1• √ ∏ ∏ × ≥ ∞≤ ≥ ∏ ∏ 1 × ∞ ≥ ≠ ∞ 2√ ≥) ≤ 1× √ √ ∏ ∏ ∞ ≤ 2≥ 1⁄ ≥ )≤ ≤ 1≥ ∏∏ ∏1 ≤ ) •≥ ° ≥ 1 ∏ ∏¬ ∏ 1 °)≠ × 1 ∏ ∏ ∏ ∏ 1 ≥)≠ × 1⁄ √ ∏ ∏ ∏ ∏ 1 √≥ °) ≥ 1 ∏ ¬ √ ∏ ¬1× ¬1 • √ )• 1≥ ∏1 ≤ εταλ ≤ ∏ 1≥ ∏ ≤ ∏2 °∏) ≠∏ ≠ ≥∏ 1 ∏ ∏ ∞ ≤ ≥ 1 × 1∞ ≤ ≥ ≠ ∞ √ ≥) 期朱建荣等 预估修正法对河口海岸海洋模式稳定性的提高海洋与湖沼 卷ΙΜΠΡΟςΕΜΕΝΤΟΦΤΗΕΣΤΑΒΙΛΙΤΨΟΦΤΗΕΕΣΤΥΑΡΙΝΕΧΟΑΣΤΑΛΟΧΕΑΝΜΟ∆ΕΛΒΨΑΠΡΕ∆ΙΧΤΟΡ−ΧΟΡΡΕΧΤΟΡΣΧΗΕΜΕ2 ≠ 2 ∏ ≥ ∏2÷ΣτατεΚεψΛαβορατορψοφΕστυαρινεανδΧοασταλΡεσεαρχη ΕαστΧηιναΝορµαλΥνιϖερσιτψ ΣηανγηαιΑβστραχτ× ∞ ∏ ≤ 2 ∞2≤ 2 ∞∏ 1× √ 2 1• 2 ≤ ∏ 2 √ √ 1 ∏ 2 ∏ ∏ ∏ ∏ ≠ ∞ ≤ ≥ 1× 2 √ √ ∏ ¬ 1× ∞ ∏ ≤ ∏ ∏ ∏ ∏ ∏ × • ≤∏ × ∏ • ≤∏ ≠ ≥ • ≤∏ ∏ 1× ∏ ∏ √ ∏ √ ∏ 1Κεψωορδσ∞ ∏ ° 2 ≤ ∏ ≠ ∞ ≤ ≥。