基于FVCOM的杭州湾三维泥沙数值模拟

Γ
= C* y， z， t） ， 本研究河 i （ x，
2
2． 1
杭州湾三维模型率定和验证
边界设定 将杭州湾 （ 上为盐官、 下为芦潮港—镇海连线 ） 剖 1 ， 20 335 个 分为如图 所示计算网格 共 12 539 个节点， 单元， 外模时间步长 1 s， 内模时间步长 3 s． 盐官和芦潮 港、 镇海均给定逐时潮位、 悬沙边界条件．
(
)
( (
) )
2 2 2 2 2 1 Dq Duq Dvq ωq q + DF + 2 D（ P + P － ε） = + + + Kq q s b D σ t x y σ σ 2 2 2 2 2 珟 Dq l Duq l Dvq l ωq l 1 W q l + DF + DlE + + + = Kq ε l 1 Ps + Pb － t x y D σ σ E1 σ
T S q l h y m m 2 2 h i 2 2 m m h h q s m b h 0 3 珟 = 1 + E 2 l2 / （ κ L2 ） ， L = （ ζ － z ） ε = q / B1 l， W －1
[
]
[
(
+ （ H + z）
－1
，
Sm =
0 ． 427 5 － 3 ． 354 G h l2 g ρ z 0 ． 494 ， Gh = ， Sh = （ 1 － 34 ． 676 G h ） （ 1 － 6 ． 127 G h ） 1 － 34 ． 676 G h k ρ0
第4 期
吴修广，等： 基于 FVCOM 的杭州湾三维泥沙数值模拟
87
wenku.baidu.com
泥沙验证． 应用模型计算了杭州湾大潮期间涨落潮泥沙输运过程 ， 通过平面、 立面的泥沙分布特征， 成功再现 ． ， 了杭州湾潮流泥沙的运动过程 尽管模拟精度达不到平面二维模型的水平 但为三维潮流泥沙模型在强潮河 口湾的应用作了重要的探索．
地形变化， 研究人员将 Phi1ips 提出的坐标变换应用到河口与海岸三维模型中
［8 ］
． 以 Princeton 大学 Mellor 为
其代表性 首的海洋动力环境数值模拟小组从 20 世纪 80 年代开始一直致力于三维数模的开发与应用研究， ECOMSED． 德国汉堡大学海洋研究所 Backhaus 等人研发的汉堡陆架海模式 HAMSOM 在世界 软件为 POM， 许多陆架海上也得到广泛应用． Sheng 建立了一般曲线坐标下的三维水动力学模型 （ CH3D ） ， 该模型也采用 S 坐标系， ． Delft Delft3D 水平方向的运动采用水平流速矢量的逆变分量来表示 荷兰 研究所建立起来 也得到 ［6 ］ ． 了较广泛应用 ECOM， FVCOM 等， 目前我国河口海岸三维数值模拟均采用国际上应用广泛的三维模型 ， 主要有 POM， ［9 ］ 其中 POM 应用时间最长、 范围最广， 几乎世界上各大海域都已应用 ． 各模型都有自身的特点和制约因素， 虽然在国、 内外不少河口海岸水域都取得过较好的计算精度 ， 但在特定海域的应用均需很长时间来探索和调 试， 我国一些学者 也 在 国 外 模 型 的 基 础 上 ， 根据我国海域特点开展了河口海岸水沙模拟及工程影响分 析
果基本满意． 应用模型计算了杭州湾大潮期间涨落潮泥沙输运过程 ， 通过平面、 立面的泥沙分布特征， 成功再现 了杭州湾潮流泥沙的运动过程 ． 尽管模拟精度达不到平面二维模型的水平 ， 但为三维潮流泥沙模型在强潮河口 湾的应用作了重要的探索 ， 通过进一步完善模型参数及分析 ， 模型可以应用于涉水建筑物对附近三维水沙影响 可为强潮河口湾开发、 水环境保护及生态建设等提供更加丰富的研究成果 ， 具有广阔的应用前景． 的分析评估，
第4 期 2011 年 12 月
水
利
水
运
工
程
学
报
HYDROSCIENCE AND ENGINEERING
No． 4 Dec． 2011
基于 FVCOM 的杭州湾三维泥沙数值模拟
吴修广，刘光生，程文龙
（ 浙江省水利河口研究院 浙江省河口海岸重点实验室 ，浙江 杭州 310020 ）
摘要： 应用 FVCOM 建立了三维潮流泥沙数学模型， 根据实测水文资料， 对杭州湾开展了三维潮流泥沙验证 ， 结
DC i DuC i DuC i C i（ ω － w C i ） 1 C i = + DF C i + + + D σ K h σ t x y σ TD TDu TDv Tω 1 T + DF + DF H ^ + + + = K T D σ h σ t x y σ SD SDu SDv Sω 1 S + DF = + + + K S D σ h σ t x y σ
［10－13 ］
．
本文应用 FVCOM 建立了杭州湾三维潮流泥沙数学模型 ， 根据实测水文资料， 对杭州湾开展了三维潮流
收稿日期： 2011－07－25 基金项目： 国家自然科学青年基金项目 （ 40806037 ） 作者简介： 吴修广（ 1974－ ） ， 男， 山东阳谷人， 高级工程师， 博士， 主要从事河口海岸潮流泥沙模拟研究 ． Email： xiuguangw@ 163． com
0 m sx
τ sy ） ， ω = 0 ，
C i = 0， σ
D S T = 0， = ［ Q （ x， y， t） － SW（ x， y， 0， t） ］ ． σ σ ρC p K h n （ 2 ） 海（ 河） 床底部（ σ = －1 ） Qb u u = D （ τ ， S = A H Dtanα （ T， ， T ， S） ； ， bx τ by ） ， ω = K h － A H S ρ0 K m Ω σ σ σ σ
88
水
利
水
运
工
程
学
报
2011 年 12 月
K m 为水平和垂向涡黏系数， K h 分别 其中： A m ， 由修正的 Mellor 和 Yamada 的 2． 5 阶湍流闭合子模型计算． A h ， ： 为水平方向和垂向热力扩散系数 A m = 0 ． 5 C0 Ω u
ζ
Ah
u
u ( 槡x ) 0． 5C Ω u = ( P 槡x )
[ ( ∫ [ ( ∫
σ
) )
] ]
( (
) )
（ 2）
Dv Duv Dv ωv + + + + Dfu = t x y σ － gD
0 1 ζ gD D v + DFy － D ρdσ' + ρσ y + D σ K m y ρ0 y σ
（ 3） （ 4） （ 5） （ 6） （ 7）
般有平面二维水流、 泥沙模型， 具有较高的精度并广泛应用
［2 ］
． 但二维模型给出的沿水深平均的潮流、 泥沙
［2－4 ］ ． 分布特征不能够完全反映钱塘江河口和杭州湾的水流 、 泥沙特点， 特别是水流、 泥沙的垂向分布 在强潮河口水域， 一方面由于地形地貌、 工程几何形态尺度以及各种障碍物的影响， 流态复杂； 另一方 面， 径流、 潮流、 波浪以及风、 柯氏力等因素是影响河口地区物质输运及沉积的主要动力条件 ， 这些动力因子 ［5－6 ］ ． 的单独或耦合作用增加了河口泥沙 、 盐分、 污染物及热量输运研究的复杂性 Leedertse 国际上对潮汐河口水动力数值模拟始于 20 世纪 60 年代后期． 基于简化过的三维浅水方程， ［7 ］ （ 1973 ） 开创性地在垂直方向采用固定分层法建立了河口 、 海湾三维潮流、 盐度模型 ． 为了更好地模拟河床
Fy ， F Ci ， FT ， FS ， Fq ， F l 的定义如下： 似壁面函数． 水平扩散项 F x ， Fx ≈
Ci
u v v + ) ] + [ 2A H ] ) ] ， F ≈ x [ 2A H( y y x y （F ， F ， F ， F， F ） ≈ [ (A H ) + (A H ) ] （ C ， T， S， q ， q l） x y x y 1 K = lqS ， K = lqS ， K = 0 ． 2 lq， P = K [ ( u ) + ( v ) ] ， P = ( gK ρ) ， ρ z z z u + 2 A H u v 2 Am H + m x y x y x
词： 杭州湾； 潮流泥沙； 三维模拟； FVCOM 中图分类号： TV148 文献标志码： A 关 键
文章编号： 1009－640X（ 2011） 04－0086－11
［1 ］ 杭州湾是典型的强潮河口湾， 潮波变形剧烈， 潮流、 泥沙运动复杂 ． 杭州湾两岸是我国经济最发达的 近年来大量的涉水工程建设改变了局部水域的潮流泥沙过程 ． 目前针对杭州湾研究的数学模型一 地区之一，
0 r
2
+
2
1 v u + 2 x y 1 v u + 2 x y
(
)
2
+ v y + v y
( )
2
+
(
)
2
( )
2
v 控制体的面积； Ω ζ 为水位 ζ 控制体的面积； P r 为 Prandtl 数． 其中： C0 为常数； Ω 为流速 u， 1 ． 2 边界条件 （ 1 ） 自由表面（ σ = 0 ） u u = D （ τ ， ) ( ρ K σ σ
(
)
(
)
C i = D i － D i ， E i = Δ tQ i（ 1 － P b ） F bi（ τ b / τ ci － 1 ） σ 式中： τ D i 为非均匀沙 i 的临界淤积切应力； τ C i 为非均匀沙 i 的临界冲刷切应力； M 为冲刷系数； 珗 τ b 为底部切应 珗 珗 ． ＜ ， ， E = w C （ 1 － / ） ； ＜ 力； 珗 τ b = ρf b V 当 τ τ 时 水中泥沙处于落淤状态 则 α τ τ 当 τ τ b ＜ τ e i 时， b Vb b Di Ci Ci i Di Di 底部泥沙边界： Kh 床面处于不冲不淤状态， 则 E C i = 0 ； 当 τ b ＞ τ e i 时， 床面泥沙起动， 则 E C i = － M（ τ / τ ei － 1 ） ． （ 3 ） 岸壁边界 Φ = 0， u， v， w， T， S， Ci ， q2 ， q2 l］ 其中 Φ = ［ n （ 4 ） 潮位、 流量、 悬沙边界河口上游边界给定潮位过程 ζ t（ 或者流量边界 ） 和温度 T t ， 盐度 S t 等的时间过 在海（ 河） 岸边界， 给定： C i （ x， 温度 T t ， 盐度 S t 等的时间过程． 入流时， 程； 外海边界给定潮位 ζ t ， y， z， t） C i C i 3 3 + u n 珗 = 0． 流边界约 0． 5 kg / m ， 湾口边界约 2 kg / m ； 出流时， t n
(
)
(
)
(
)
（ 8） （ 9）
S） ρ = ρ（ T，
D 分别为潮位和动态水深； H 为静态水深； u， v， y， 式中： σ 为垂向相对坐标； t 为时间； ζ， ω 分别为 x， σ 方向上 T 分别为盐度和温度； g 为重力加速度； f 为柯氏力参数（ f = 2 Φsinφ， 的流速分量； C i 为非均匀泥沙 i 的浓度； S，
2 l 分别为湍流动能和湍流宏观尺度 ； K m ， Kh φ 为纬度， Φ 为地球自转角速度） ； ρ， ρ0 分别为海水和淡水密度； q ， ^ 分别为垂向涡黏系数和垂向热力扩散系数 ； K q 为湍流动能垂向涡黏扩散系数； w C i 为非均匀沙 i 的沉速； H 为 珟为近 P b 分别为剪切和浮力引起的湍流动能产生项 ； ε 为湍流动能耗散率； W 太阳短波辐射的垂向梯度项； P s ，
1
1． 1
数学模型
模型控制方程 ζ Du Dv ω + + + =0 t u y σ
2 Du Du Duv uω + + + － Dfv = t x y σ
（ 1）
－ gD
0 1 ζ gD D u + DFx － D ρdσ' + ρσ x + D σ K m x ρ0 x σ σ 2