Page 55 - 2023年第54卷第3期

图 １ 网格划分及插值示意

向耦合。首先设置水深阈值，根据该阈值将细网格分为淹没区和非淹没区。若网格内的水深小于阈
值，将其视为非淹没区，采用地表水文模型计算；反之，视为淹没区，采用二维水动力模型计算。随
着降雨条件的变化，网格水深不断发生变化，淹没区和非淹没区范围也发生变化，采用 ＣＭＩ连接淹没
区和非淹没区。
动态双向耦合模型以 Ｇｏｄｕｎｏｖ型通量求解分析为理论基础，根据接触间断处的水位、地表高程及
特征波波速计算界面通量及两种模型之间的水量、动量交换。图 ２展示了 ｘ方向的耦合动边界的移动
情况，在 ｘ方向耦合界面左右两侧的特征波如下：
Ｓ ＝ ｕ － ｇｈ ｉ，ｊ （４）
槡
Ｌ
ｉ，ｊ
Ｓ ＝ ｕ ｉ ＋ １，ｊ ＋ ｇｈ ｉ ＋ １，ｊ （５）
槡
Ｒ
以下根据特征波的大小，对耦合边界的不同情况详细说明。
第一种情况，两个特征波均大于 ０，即 Ｓ≥Ｓ＞０，耦合边界位置向淹没区移动，如图 ２（ａ）所示。
Ｌ
Ｒ
此时耦合界面的通量计算仅依赖于网格 （ ｉ，ｊ），与网格（ｉ ＋ １ ，ｊ）无关。采用地表水文模型计算网格
（ ｉ，ｊ）的水深和流速，采用水动力模型计算网格（ｉ ＋ １ ，ｊ）的水深和流速，在耦合界面上只有水量通过。
网格（ ｉ ＋ １，ｊ）的速度不影响网格（ｉ，ｊ）。通过耦合界面的通量计算如下：
Ｆ ＝ Ｆ ＝ ［ｈｕ，０］ Ｔ （６）
ｉ ＋ １?２，ｊ Ｌ ｉ，ｊ
从上式可以看出，此时在耦合界面上只传递质量，不传递动量。
第二种情况，两个特征波一个大于 ０，一个小于 ０，即 Ｓ＜０＜Ｓ，此时耦合界面位置不移动，如图
Ｌ Ｒ
２（ｂ）所示。耦合界面的通量计算与网格（ｉ，ｊ）和网格（ｉ ＋ １，ｊ）均有关，计算公式如下：
ＳＦ － ＳＦ ＋ ＳＳ（Ｕ － Ｕ ）
Ｌ
Ｒ
Ｌ
Ｌ Ｒ
Ｌ
Ｒ
Ｒ
Ｆ ＝ ｆ（Ｆ，Ｆ ） ＝ （７）
ｉ ＋ １?２，ｊ Ｌ Ｒ
Ｓ－ Ｓ Ｌ
Ｒ
Ｔ
２
２
Ｔ
式中：Ｕ ＝ ［ｈ，０］ ；Ｕ ＝ ［ｈ，ｈｕ］ Ｔ ；Ｆ ＝ ［ｈｕ，０］ ；Ｆ ＝ ［ｈｕ，ｈｕ＋ ｇｈ?２］ Ｔ 。
ｉ ＋ １ ，ｊ
ｉ，ｊ
Ｒ
ｉ ＋ １ ，ｊ
Ｌ
Ｒ
Ｌ
ｉ，ｊ
第三种情况，两个特征波均小于 ０，即 Ｓ＜０，Ｓ＜０，此时耦合界面位置向非淹没区移动，如图 ２（ｃ）
Ｒ
Ｌ
所示。耦合界面的通量计算仅与网格（ ｉ ＋ １，ｊ）有关，计算如式（８）所示。
２
２
Ｆ ｉ ＋ １?２，ｊ ＝ Ｆ ＝ ［ｈｕ，ｈｕ＋ ｇｈ?２］ Ｔ （８）
ｉ ＋ １，ｊ
Ｒ
从式（７）（８）可以看出，在这两种情况下，网格（ｉ ＋ １，ｊ）的流速对网格（ｉ，ｊ）产生影响，耦合界面
处既有质量传递又有动量传递，通过耦合边界的质量和动量守恒。
同理，在 ｙ方向上，耦合界面的变化亦可采用相同的计算方法。
在 ｙ方向耦合界面上下两侧的特征波为：
Ｓ ＝ ｖ ＋ ｇｈ （９）
ＵＰ ｉ，ｊ ＋ １ 槡 ｉ，ｊ ＋ １
Ｓ ＝ ｖ － ｇｈ （１０）
ＤＮ ｉ，ｊ 槡 ｉ，ｊ
— ３ ０ ５ —