前言
CFD小工坊是一个本博客一个全新的系列，他将包含我学习二维浅水方程理论，直至编译一个实用的求解器的过程。主要内容有：
- 二维浅水方程基本方程，物理意义；
- 数值离散方法（非结构化网格），以及离散方程的形式；
- 求解器的编写。
对于编程语言，我准备以c/c++为主，并辅以一定的MATLAB脚本。我目前还未熟练掌握c/c++，也希望通过这个系列能逐步学习其使用。也希望和大家一起交流，不断完善我们的博客。
对于本博客中所涉及的所有源代码，我目前暂不开源，但我会在博客中展示一些核心部分。待整个模型开发完成且手册编写完成后，我大概会将源代码和相关资料上传至GitHub。
作为这个系列的开始，本篇博文将简要介绍浅水方程。
浅水动力学
顾名思义，浅水动力学的研究对象是浅水流。浅水流一般有着以下特征：
- 有自由表面；
- 以重力为主要驱动力，以水流与面体边界之间及水流内部的摩阻力为主要耗散力，有时还存在水面气压场、风应力及地转偏向力的作用；
- 水平流速沿垂线近似均匀分布，不必考虑实际存在的对数或指数等形式的垂线流速分布;
- 水平运动尺度远大于垂直运动尺度，垂向流速及垂向加速度可忽略，垂向的水压力分布接近静压分布。
事实上，自然界中的水流运动在物理上都属于三维流动，其运动一般可用不可压缩Navier-Stokes方程描述。考虑到浅水流动有着上述第3、4点特性，我们通常将其动力学方程沿着垂向积分（忽略物理量垂向上的空间变化），得到一个二维方程组。这个描述浅水流动的方程组被称为浅水方程，或圣维南方程组（2-D Saint-Venant Equation, 2D SVE）。
在自然界和工程应用中，能作为浅水流动处理的水流通常出现在下列情况中：
- 水深相对较浅，即流动的水平尺度远大于垂直尺度。
- 水底底坡较缓；如果底坡倾角为α，则有α≈sinα≈tanα。此时，底坡引起的垂直速度和垂直环流可以忽略，也不必考虑垂直加速度及由此产生的动水压力（满足“静压假定”）。
- 水面渐变且坡度较缓。
- 无明显的垂直环流。
因此，浅水方程常用于描述以下流动：漫滩洪水、灌溉水系、平原河网水流、入库洪水流、溃坝流、河口潮流、近岸风暴潮等等。
然而，从数学上求解二维浅水方程并不容易。通常需要数值方法将其离散，求其近似解。这边衍生出“计算浅水动力学”这一学科方向。至今，已有许多求解浅水方程的成熟方法，本博客将选择其中的一种“有限体积法”来设计、完成求解器。
二维浅水方程
二维浅水方程的守恒形式如下：
∂
U
∂
t
+
∂
E
(
U
)
∂
x
+
∂
G
(
U
)
∂
y
=
S
(
U
)
U
=
(
h
h
u
h
v
)
,
E
(
U
)
=
(
h
u
h
u
2
+
g
h
2
2
h
u
v
)
,
G
(
U
)
=
(
h
v
h
u
v
h
v
2
+
g
h
2
2
)
,
S
(
U
)
=
(
0
g
h
(
S
0
x
−
S
f
x
)
g
h
(
S
0
y
−
S
f
y
)
)
\dfrac{\partial \bold{U}}{\partial t} + \dfrac{\partial \bold{E(U)}}{\partial x} + \dfrac{\partial \bold{G(U)}}{\partial y} = \bold{S(U)} \\[6pt] \bold{U} = \left( \begin{matrix} h \\ hu \\ hv \end{matrix} \right), \bold{E(U)} = \left( \begin{matrix} hu \\ hu^2+\dfrac{gh^2}{2} \\ huv \end{matrix} \right), \bold{G(U)} = \left( \begin{matrix} hv \\ huv \\ hv^2+\dfrac{gh^2}{2} \end{matrix} \right), \\[6pt] \bold{S(U)} = \left( \begin{matrix} 0 \\ gh(S_{0x} - S_{fx}) \\ gh(S_{0y} - S_{fy}) \end{matrix} \right)
∂t∂U+∂x∂E(U)+∂y∂G(U)=S(U)U=
hhuhv
,E(U)=
huhu2+2gh2huv
,G(U)=
hvhuvhv2+2gh2
,S(U)=
0gh(S0x−Sfx)gh(S0y−Sfy)
式中，x和y表示平面笛卡尔坐标的两个方向，t表示时间，h表示水深；u和v表示x和y方向上的流速，**E(U)和G(U)**表示界面通量，**S(U)**表示源项。
S
0
x
=
−
∂
z
b
∂
x
,
S
0
y
=
−
∂
z
b
∂
y
S
f
x
=
n
2
u
u
2
+
v
2
h
−
4
/
3
,
S
f
y
=
n
2
v
u
2
+
v
2
h
−
4
/
3
S_{0x} = -\dfrac{\partial z_b}{\partial x}, S_{0y} = -\dfrac{\partial z_b}{\partial y} \\[6pt] S_{fx} = n^2 u \sqrt{u^2+v^2} h^{-4/3}, S_{fy} = n^2 v \sqrt{u^2+v^2} h^{-4/3}
S0x=−∂x∂zb,S0y=−∂y∂zbSfx=n2uu2+v2h−4/3,Sfy=n2vu2+v2h−4/3
式中，*zb*表示水底高程，n表示水底糙率。
此外，需要明确的一点是水位η、水深h及水底高程zb的相对关系是：
η
=
h
+
z
b
\eta = h + z_b
η=h+zb
参考资料
- 谭维炎《计算浅水动力学》
- Jingming Hou et. al. A 2D well-balanced shallow flow model for unstructured grids with novel
slope source term treatment. Advances in Water Resources. 2013.
