浅水波方程是假设垂直方向尺度远远小于水平方向尺度时从三维Navier-Stokes方程推导而出的二维或一维空间下的流体运动方程。广泛地应用与海啸,大气环流,河道输运,虚拟现实等问题。
如此我们的目标是要将Navier-Stokes方程
应用到图1所示问题。它有一个底面b和流体的上表面 ζ 。流体的高度为H,它是b和 ζ 的和。
注意到方程(1)中的b为流体所受体积力。在浅水波方程中最常见的体积力为科氏力 (k为垂直方向矢量)。这一项在研究洋流,大气环流等中是必不可少的。但是为简便起见,下文中我们将不考虑体积力。
在如上情况下,静水压力p是可以直接求出的。设流体上表面的压力为 ps ,则
将其带回方程式(1),得
下面对方程式(2)沿高度方向作平均计算。取流速u的平均值为
(注:对z向速度的处理方法不同会得到不同的方程【1】,另外一种常见的浅水波方程是Boussinesq方程,该方程考虑了沿z向的速度分布),则(1)式中的连续方程化简为
动量平衡方程(2)化简为
上式中的f为流体上下表面(空气和河床)所受粘性力
方程(3),(4)为本次待求的浅水波方程方程。
为了方便有限元求解,将方程(3),(4)改写为
其弱形式为
这就是我们要使用的有限元控制方程。
考虑2000m长的水库在中间被大坝隔断。水库上游初始高度10m,下游高度5m。假设水坝在一瞬消失。该浅水波流动问题有理论解【2】,常被由于检验软件计算的正误和算法的精度。
图2和图3分别为采用开源软件(软件名待定)的计算结果和参考文献【3】=理论解的计算结果。这个参考计算结果实在是太好,跟理论解几乎一致。本计算尝试了RK Backward Euler,DIRK 1 Stage Theta Method, RK Implicit 2 Stage 2nd order Lobatto IIIB, BDF2等时间积分格式都未得到理想的效果。可见常套算法对此方程不太好使。
溃坝过程模拟结果
【1】
https://arxiv.org/pdf/1706.08815.pdfarxiv.org/pdf/1706.08815.pdf
【2】 J.J. Stoker. Water waves: the Mathematical Theory with Applications. Wiley Interscience, 1992
【3】
http://etd.lib.metu.edu.tr/upload/12618378/index.pdfetd.lib.metu.edu.tr/upload/12618378/index.pdf