一种浸入边界法的壁面函数湍流计算方法与流程

1.本发明涉及计算流体力学数值模拟技术，尤其是一种基于浸入边界法的壁面函数湍流计算方法。


背景技术：

2.在自然界或者工程领域，很多问题涉及物体大变形下的流固耦合问题，比如鸟类或昆虫挥动翅膀飞行，降落伞或翼伞的展开过程等。对于这类问题，一方面物体经历了很大的结构变形，引起周围流场的变化；另一方面，周围流场的变化又会反过来影响物体的变形，因此物体经历了复杂的流固耦合过程。流固耦合问题的模拟一直以来都是计算科学领域的经典问题，也是一个难点。传统贴体网格方法能够较好地表示物体表面，但是在处理流固耦合问题时，需要采用动网格方法改变流场网格，以保证网格质量。当物体经历大变形时，流场网格质量往往会急剧下降，甚至出现负体积，导致计算失败。浸入边界法是另一种流场求解方法，其特点是空间网格固定，即使物体发生了巨大变形，空间网格也固定不动，因此避免了网格质量下降甚至负体积的问题，是处理大变形流固耦合问题的一种有效方法，同时方便自动生成网格。但是，浸入边界法的代价是壁面边界处理较贴体网格方法更为复杂，需要采用特殊方法，并且边界处的积分量计算也较贴体网格方法更为复杂。
3.自然界和工程中的流动问题一般是高雷诺数问题。对于这类问题，传统贴体网格方法一般采用长宽比极大的网格单元，以便准确捕捉边界层效果。但是对于浸入边界法，一般采用各向同性的笛卡尔网格进行计算。因此为了准确捕捉边界层，在物面附近所需的网格数目巨大，极大地影响了浸入边界法在高雷诺数流动问题的应用。


技术实现要素：

4.为了解决高雷诺数下浸入边界法计算湍流时物面附近区域所需网格数巨大的问题，本发明提供了一种基于spalart-allmaras壁面函数的可压缩流浸入边界法湍流计算方法，旨在通过壁面函数减少湍流计算时物面附近区域所需的网格数，提高浸入边界法处理湍流流动的能力。
5.为了实现上述目的，本发明提供的技术方案是：一种浸入边界法的壁面函数湍流计算方法，包括：
6.步骤1：建立主流控制方程，所述主流控制方程为favre-averaged navier-stokes方程，
7.步骤2：建立湍流模型方程，所述湍流模型为spalart-allmaras湍流模型方程；步骤3：根据所述主流控制方程，在远离壁面的区域采用格点格式有限体积法进行空间离散，而在壁面附近区域采用fiver浸入边界法进行处理，建立如下的一维半riemann问题：
8.9.其中，为物面附近的流场网格点沿着物面局部法向方向的流场状态,ρ是密度，ξ是一维局部法向方向坐标，指向物体内部，即沿着方向，是物体运动速度，δtn表示当前时间步；
10.步骤4：根据所述一维半riemann问题的精确解计算壁面处的流场状态，从而得到壁面处的通量；
11.步骤5：根据所述湍流模型方法，在远离壁面的区域进行空间离散，在壁面附近区域采用fiver浸入边界法计算密度、法向速度和总能，特别采用改进的壁面函数计算切向速度和涡粘性系数。
12.作为优选的一种技术方案，所述步骤5进一步包括：
13.步骤51：由壁面附近通量计算位置和最近壁面的法向确定流场内部点；
14.步骤52：以该流场内部点到壁面的距离为截断距离，修改spalart-allmaras湍流模型方程中的衰减函数；
15.步骤53：将该点相对壁面的切向速度值代入壁面函数，通过确定壁面的摩擦速度；
16.步骤54：将壁面的摩擦速度和通量计算位置的壁面距离代入壁面函数得到其切向速度，由壁面的摩擦速度和壁面距离得到涡粘性系数；
17.步骤55：计算得到通量计算位置处的通量。
18.作为优选的一种技术方案，所述流场内部点的流场状态值由其所在单元通过线性插值得到。
19.作为优选的一种技术方案，所述流场内部点到壁面的距离为壁面附近网格尺寸的3倍。
20.作为优选的一种技术方案，所述步骤3中的空间网格采用的是四面体单元，采用格点格式有限体积法离散，流场状态存储在网格单元的节点上，且空间网格不含有物体表面网格，仅仅在空间网格生成时在壁面附近区域进行局部网格加密以便准确捕捉物体附近流动状态；
21.作为优选的一种技术方案，在所述湍流模型方程中，在远离物面区域采用一阶迎风格式计算其对流项，采用中心差分计算其耗散项，在靠近壁面区域，截断距离以内spalart-allmaras湍流模型方程采用修改的衰减函数。
22.作为优选的一种技术方案，所述步骤s53中的摩擦速度由牛顿迭代法求解。本发明相对于现有技术的有益效果是：
23.(1)本发明通过结合spalart-allmaras壁面函数和fiver浸入边界法实现高雷诺数湍流问题的浸入边界法数值模拟计算，减少了湍流计算所需的边界层网格数；
24.(2)本发明基于可压缩navier-stokes方程，可处理高速流动问题；
25.(3)本发明对于物体形状以及是否封闭没有限制，具有广泛的适用性。
附图说明
26.图1是本发明一实施例提供的壁面附近区域示意图。
具体实施方式
27.以下结合附图对本发明的实施例作进一步详细描述。
28.本实施例提供一种基于spalart-allmaras壁面函数的可压缩流浸入边界法湍流计算方法，减少了浸入边界法在高雷诺数下计算湍流问题所需的壁面附近网格数量，提高了浸入边界法计算高雷诺数下湍流问题的能力。
29.为了解决上述问题，本实施例提供一种浸入边界法的壁面函数湍流计算方法，具体包括一下步骤：
30.步骤1：建立主流控制方程，在本实施例中，主流控制方程为favre-averagednavier-stokes方程，主流控制方程表示为
[0031][0032]
其中，w为守恒变量(密度，动量和总能)，ω为控制体的体积(三维)或面积(二维)，为控制体边界，fc表示对流通量，fv表示粘性通量。当考虑湍流问题时，表示为favre-averaged navier-stokes方程，流场变量为平均值。
[0033]
步骤2：建立湍流模型方程，在本实施例中，湍流模型方程采用spalart-allmaras湍流模型方程。
[0034]
步骤3：对于主流控制方程，在远离壁面的区域采用格点格式有限体积法进行空间离散，对流通量计算采用roe格式，由线性重构实现二阶精度，粘性通量计算采用中心差分，其中流场变量的梯度由green-gauss方法获得。壁面附近区域采用fiver(finite volume method with exact two-material riemann problems)浸入边界法处理壁面边界条件，建立如下的一维半riemann问题：
[0035][0036]
其中，为物面附近的流场网格点沿着物面局部法向方向的流场状态,ρ是密度，ξ是一维局部法向方向坐标，指向物体内部，即沿着方向，是物体运动速度，δtn表示当前时间步。
[0037]
步骤4：根据所述一维半riemann问题的精确解计算壁面处的流场状态，从而得到壁面处的通量。
[0038]
在本实施例中，该一维半riemann问题的精确解分为两种情况，一种是激波，另一种是稀疏波。再由精确解计算壁面处的流场状态(密度、法向速度和总能)，对于无粘或者层流问题，由roe格式计算得到壁面处的通量。如图1中曲线表示物体边界，i点为流场内部点，j点为物体内部点，在边ij的中点c处求解一维半riemann问题，得到点c处的密度，法向速度和总能。
[0039]
步骤5：对于spalart-allmaras湍流模型方程，在远离壁面区域采用类似主流方程的方法离散spalart-allmaras湍流模型方程，在壁面附近区域采用fiver浸入边界法获得密度、法向速度和总能，采用改进的壁面函数计算壁面附近区域的切向速度和涡粘性系数，其步骤为：
[0040]
a.由壁面附近通量计算位置和最近壁面的法向确定流场内部一个点，该流场内部点的流场状态值由其所在单元通过线性插值得到，如图1所示p点，由c点沿着最近壁面的法向方向得到，p点到壁面的距离为壁面附近网格尺寸的3倍。p点的流场状态由p点所在的单元通过线性插值计算。
[0041]
b.以该流场内部点到壁面的距离为截断距离，修改spalart-allmaras湍流模型方程中的衰减函数f
v1
，即当流场内部点到壁面的距离小于p点到壁面的距离时，衰减函数为
[0042][0043]
其中，rd＝d
p
/d，d
p
为p点到壁面的距离，d为流场内部点到壁面的距离，c
v1
＝0.71，c为涡性系数，ν为运动学粘性系数。
[0044]
c.spalart-allmaras壁面函数为
[0045][0046]
在本实施例中，a1＝8.14822158,a2＝-6.92870938,b1＝7.46008761,b2＝7.46814579,c1＝2.54967735,c2＝1.33016516,c3＝
[0047]
3.59945911,c4＝3.63975319，u
+
＝u/u
τ
，y
+
＝u
τ
d/ν，u
τ
为壁面摩擦速度。
[0048]
对于切向速度，有如下关系
[0049]ut
＝u
τ
f(y
+
)(5)
[0050]
函数f表示式(3)。将p点的切向速度式(5)，通过牛顿迭代法确定壁面摩擦速度u
τ
；
[0051]
d.将u
τ
和c点的壁面距离代入式(5)得到c点的切向速度和涡粘性系数
[0052][0053]
e.以c点处的流场状态为右状态，计算c点的通量，从而实现主流方程和spalart-allmaras湍流模型方程的壁面边界条件。
[0054]
本发明相对于现有技术的有益效果是：
[0055]
(1)本发明通过结合spalart-allmaras壁面函数和fiver浸入边界法实现高雷诺数湍流问题的浸入边界法数值模拟计算，减少了湍流计算所需的边界层网格数；(2)本发明基于可压缩navier-stokes方程，可处理高速流动问题；
[0056]
(3)本发明对于物体形状以及是否封闭没有限制，具有广泛的适用性。
[0057]
以上仅是本发明的优选实施方式，本发明的保护范围并不仅局限于上述实施例，凡属于本发明思路下的技术方案均属于本发明的保护范围。应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理前提下的若干改进和润饰，应视为本发明的保护范围。

技术特征：
1.一种浸入边界法的壁面函数湍流计算方法，其特征在于，包括：步骤1：建立主流控制方程，所述主流控制方程为favre-averaged navier-stokes方程，步骤2：建立湍流模型方程，所述湍流模型为spalart-allmaras湍流模型方程；步骤3：根据所述主流控制方程，在远离壁面的区域采用格点格式有限体积法进行空间离散，而在壁面附近区域采用fiver浸入边界法进行处理，建立如下的一维半riemann问题：其中，为物面附近的流场网格点沿着物面局部法向方向的流场状态,ρ是密度，ξ是一维局部法向方向坐标，指向物体内部，即沿着方向，是物体运动速度，δt
n
表示当前时间步；步骤4：根据所述一维半riemann问题的精确解计算壁面处的流场状态，从而得到壁面处的通量；步骤5：根据所述湍流模型方法，在远离壁面的区域进行空间离散，在壁面附近区域采用fiver浸入边界法计算密度、法向速度和总能，特别采用改进的壁面函数计算切向速度和涡粘性系数。2.根据权利要1所述的壁面函数湍流计算方法，其特征在于，所述步骤5进一步包括：步骤51：由壁面附近通量计算位置和最近壁面的法向确定流场内部点；步骤52：以该流场内部点到壁面的距离为截断距离，修改spalart-allmaras湍流模型方程中的衰减函数；步骤53：将该点相对壁面的切向速度值代入壁面函数，通过确定壁面的摩擦速度；步骤54：将壁面的摩擦速度和通量计算位置的壁面距离代入壁面函数得到其切向速度，由壁面的摩擦速度和壁面距离得到涡粘性系数；步骤55：计算得到通量计算位置处的通量。3.根据权利要求1所述的壁面函数湍流计算方法，其特征在于，所述流场内部点的流场状态值由其所在单元通过线性插值得到。4.根据权利要求2或3所述的壁面函数湍流计算方法，其特征在于，所述流场内部点到壁面的距离为壁面附近网格尺寸的3倍。5.根据权利要求1所述的壁面函数湍流计算方法，其特征在于，所述步骤3中的空间网格采用的是四面体单元，采用格点格式有限体积法离散，流场状态存储在网格单元的节点上，且空间网格不含有物体表面网格，仅仅在空间网格生成时在壁面附近区域进行局部网格加密以便准确捕捉物体附近流动状态。6.根据权利要求1所述的壁面函数湍流计算方法，其特征在于，在所述湍流模型方程中，在远离物面区域采用一阶迎风格式计算其对流项，采用中心差分计算其耗散项，在靠近壁面区域，截断距离以内spalart-allmaras湍流模型方程采用修改的衰减函数。7.根据权利要求2所述的壁面函数湍流计算方法，其特征在于，所述步骤s53中的摩擦
速度由牛顿迭代法求解。

技术总结
本发明公开了一种浸入边界法的壁面函数湍流计算方法，涉及计算流体力学领域，能够提高FIVER浸入边界法求解高雷诺数湍流问题的计算效率。本发明包括：结合Spalart-Allmaras壁面函数和FIVER浸入边界法建立一种基于壁面函数的湍流计算方法，利用壁面函数计算壁面附近网格的切向速度和修改的涡粘性，再确定壁面附近通量计算所需的相关变量。本发明减少了FIVER浸入边界法在高雷诺数下计算湍流问题所需的壁面附近网格数量，提高了FIVER浸入边界法计算高雷诺数下湍流问题的能力。法计算高雷诺数下湍流问题的能力。法计算高雷诺数下湍流问题的能力。


技术研发人员：高宜胜 张红英 李正达 陆伟伟 刘琦 张芃 童明波
受保护的技术使用者：航宇救生装备有限公司
技术研发日：2023.03.21
技术公布日：2023/7/17