[发明专利]一种考虑非定常摩擦模型的自由表面瞬态流的FVM模拟方法

权利要求书

1.一种考虑非定常摩擦模型的自由表面瞬态流的FVM模拟方法，其特征在于，包括如下步骤：

S1：构建包含布鲁诺非定常摩擦模型的管道自由表面瞬态流动控制方程；

S2：根据有限体积法的计算特性，将计算区域进行网格划分；

S3：基于步骤S2的网格划分，采用有限体积法计算包含了布鲁诺非定常摩擦模型的自由表面流动；

S4：根据步骤S3的计算结果，进行模型评价验证。

2.根据权利要求1所述的一种考虑非定常摩擦模型的自由表面瞬态流的FVM模拟方法，其特征在于，所述步骤S1中引入布鲁诺非定常摩擦模型的管道自由表面瞬态流动控制方程为：

式中变量W，通量G，和源项S可以写成:

式中，Area为水的截面积；Q是流量；为截面上水柱的平均压力；ρ是液体的密度；g是重力加速度；Slope为管道的斜率；S_f是管道摩擦力；

将改进的布鲁诺非定常摩擦模型引入控制方程，具体模型如下：

S_f＝S_q+S_b#(3)

式中，S_q为稳态摩擦力项；S_b是非定常摩擦项。

3.根据权利要求2所述的一种考虑非定常摩擦模型的自由表面瞬态流的FVM模拟方法，其特征在于，所述步骤S1中根据曼宁公式，稳态摩擦力项S_q可定义为：

式中，n为管道曼宁粗糙度系数，R为水力半径；

改进的鲁诺非定常摩擦模型提出了以下计算非定常摩擦项S_b的公式：

式中，c为波速；若v0，则SGN(v)＝-1，若v≥0，SGN(v)＝1；f为摩阻系数。

4.根据权利要求3所述的一种考虑非定常摩擦模型的自由表面瞬态流的FVM模拟方法，其特征在于，所述步骤S1中利用瓦迪剪切衰减系数B*解析地预测f：

5.根据权利要求1所述的一种考虑非定常摩擦模型的自由表面瞬态流的FVM模拟方法，其特征在于，所述步骤S2中将计算区域进行网格划分的具体方式为：将整个计算域进行有限体积网格划分，建立管道自由表面流动有限体积体系的计算网格，将自由表面流管道空间域离散为N个单元，每个单元长度为Δx，时间域离散为时间间隔Δt，控制体中点编号为i，第i个控制体左右边界分别为i-1/2和i+1/2，计算模型局部离散网格。

6.根据权利要求1所述的一种考虑非定常摩擦模型的自由表面瞬态流的FVM模拟方法，其特征在于，所述步骤S3中有限体积法求解时的积分具体为：

将式(1)沿空间域x方向从边界i-1/2到边界i+1/2进行积分，并带入得到：

式中：上标n+1为t+Δt时刻；n为t时刻。

7.根据权利要求6所述的一种考虑非定常摩擦模型的自由表面瞬态流的FVM模拟方法，其特征在于，所述步骤S3中包含了布鲁诺非定常摩擦模型的自由表面流动的计算具体包括如下步骤：

A1：二阶精度重构

采用中心迎风格式方法进行二阶线性重构，从而获得空间为二阶精度FVM格式，同时引入MINMOD斜率限制器函数，避免线性插值重构时出现虚假振荡的现象；

二阶精度重构具体包括如下步骤：

A1-1：重构

式中：Ψ_jⁿ＝(W_j+1ⁿ-W_jⁿ)/Δx，Ψ_j-1ⁿ＝(W_jⁿ-W_j-1ⁿ)/Δx；

A1-2：时间推进

A2：双激波假设求解

应用双激波假求解器进行近似求解，结合守恒定律，给出了双激波假求解器中的通量G_i+1/2；

式中，O_L、O_R分别为左右测的波速；G_L和G_R分别为计算单元交界面左侧和右侧的通量。

W_L和W_R是界面左右测的变量；

O_L＝v_L-E_L#(11)

O_R＝v_R+E_R#(12)

式中，E_K(K＝L，R)可表示为：

其中Area_*为星区面积的估计值；为计算星区过水断面形心点压强；为计算星区左右两侧过水断面形心点压强；c_K(K＝L，R)为计算星区左右两侧重力波的速度；对于圆形通道，变量c、Area和可以表示为：

其中，d为管道直径；φ为水面的圆心角；h是水深；

对于Area_*的求解，基于以下公式对Area_*的值进行估计：

Area_*＝0.5[1+(V_L-V_R)/(γ_L+γ_R)](Area_L+Area_R)#(17)

其中γ＝∫(c/A)dA；对于圆截面管道，γ可由如下表示：

A3：边界条件求解

所有计算域均采用式(7)计算，包括上下游边界G_1/2和G_N+1/2；对于每个G_i+1/2，使用中心迎风格式方法进行线性重建需要接口左右两侧的数据，在计算域两侧分别加入2个计算单元，即总共加入4个计算单元，这四个计算单元的所有参数值都与原始边界处的值相同；

基于边界处理方法计算域由1～N，拓展为-1～N+2，前一时刻计算域-1～N+2的所有变量都已知，对计算域1～N的元素进行重构，可以得到上下边界的变量值；

对于圆截面通道，如果在上游边界处规定入流，和为与有关的未知变量；因此，式(19)中只有一个未知变量通过迭代直接求出和下游边界也可以用同样的方法找到；

为了提高边界的精度，边界两侧引入了2个计算单元，其参数与边界参数一致，使在提高精度的同时不引入其他误差，给出了边界单元格具体表达式：

A4：考虑布鲁诺非定常摩擦模型的时间积分

G_i-1/2和G_i+1/2的二阶格式在前三步中已求得，为了将解从时间n推进到时间n+1，需要进行时间积分，在存在布鲁诺非定常摩擦的情况下，进行以下时间积分，将布鲁诺非定常摩擦模型引入求解

A5：稳定性约束

稳定性约束包括库朗数约束，考虑了布鲁诺非定常摩擦模型的源项约束

整体的稳定性约束为：

Δt_max＝min_i＝1....n[Δt_S,max,Δt_Cr,max]#(28)。