一种内孤立波影响下潜艇姿态预测方法

本发明涉及一种内孤立波影响下潜艇姿态预测方法，属于水下航行器安全航行研究领域。

背景技术：

1、内孤立波是层状海洋中普遍存在的现象，具有中细尺度、振幅大、流速强等特征。内孤立波产生的剪切流和湍流有助于全球能量在海洋中的传输和耗散，但是也对生物活动、近海工程和潜艇航行等产生影响。由于内孤立波在海面以下发展，传播频率低、距离远，使得其观测困难，严重威胁海洋工程结构和水下航行器的安全。上世纪，美国“长尾鲨”号核潜艇，以色列“达卡”号潜艇都遭受内波影响而发生沉没。而2021年，印尼“南伽拉”号同样受到内波影响而发生沉没。由此可见，悄无声息的内孤立波已经严重影响潜艇安全航行。

2、针对上述问题，国内外均开展了相关研究，目前主要采用试验方法和数值方法开展内孤立波影响下水下结构姿态变化。前者主要依赖水池试验，采用水平推板、双推板、速度入口等方法开展研究，但上述试验方法均采用缩比模型开展研究，受试验硬件影响，模拟不同深度和不同幅度的内孤立波较为困难。试验手段无法开展波长、幅度达数百米的内孤立波模拟，同时，对缩比模型自由运动姿态数据采集也较为困难。而数值方法，同样采用与试验相同原理制造内孤立波，虽然可以开展等比潜艇模型开展研究，但波长长达数百米，持续时间达数十分钟，并且无法考虑海床地貌对内孤立波的影响，同时对计算硬件提出了更高要求，非常不利于开展内孤立波影响下潜艇姿态研究。

3、综上所述，现有试验方法和数值方法都不利于开展内孤立波影响下潜艇姿态研究，为后续潜艇安全航行提出对应解决方案。因此，水下航行器安全航行研究领域急需一种内孤立波影响下潜艇姿态预测方法。

技术实现思路

1、针对上述问题，本发明所要解决的技术问题是提供一种内孤立波影响下潜艇姿态预测方法。

2、本发明为解决上述技术问题采用以下技术方案：

3、一种内孤立波影响下潜艇姿态预测方法，所述预测方法包括以下具体步骤：

4、步骤一：建立流场域几何模型，所述流场域几何模型长度方向为内孤立波传播方向，高度方向为海洋深度方向，同时对所述流场域几何模型进行网格划分，建立流场域网格模型；

5、步骤二：定义流场域网格模型局部坐标系xfluid_local，所述流场域局部坐标系xfluid_local的原点为流场域网格模型形心，x、y和z轴方向分别为流场域网格模型长度方向、高度方向和宽度方向；

6、步骤三：定义位于内孤立波传播方向下游的流场域网格模型端面和流场域网格模型上下面为速度入口型边界条件，定义位于内孤立波传播方向上游的流场域网格模型端面为压力出口型边界条件，定义流场域网格模型剩余两面为对称型边界条件；

7、步骤四：建立潜艇及围绕潜艇的随动流场域几何模型，同时对所述潜艇及围绕潜艇的随动流场域几何模型进行网格划分，建立潜艇及围绕潜艇随动流场域网格模型；

8、步骤五：定义潜艇船体坐标系xsubmarine_local，所述潜艇船体坐标系xsubmarine_local的原点为潜艇质心，x、y和z轴方向分别为潜艇长度方向、高度方向和侧向；

9、步骤六：定义潜艇及围绕潜艇随动流场域网格模型中的潜艇表面为壁面型边界条件；

10、步骤七：根据潜艇初始姿态和潜艇相对内孤立波初始位置，将潜艇及围绕潜艇随动流场域网格模型嵌入流场域网格模型，建立总体装配网格模型；

11、步骤八：对所述总体装配网格模型进行初始化；

12、步骤九：在第n个计算时间步，对局部坐标系xfluid_local下的流场域网格模型进行整体位移，其中，整体位移矢量dfluid在局部坐标系xfluid_local下的x、y和z轴方向的三个分量分别为c·δt、0和vsubmarine_z·δt，c为内孤立波波速，vsubmarine_z为海域全局坐标系x下潜艇质心运动速度矢量v在局部坐标系xfluid_local中z向投影分量，δt为计算时间步长；这里需要说明的是，海域全局坐标系x是以整个海域为参考，以任意一地理坐标(xx经度，xx纬度)为原点，地理东西方向为全局坐标系x的x轴，地理南北方向为z轴，深度方向为y轴(上述步骤七所述潜艇初始姿态，比如艇艏朝向就是以海域全局坐标系x进行描述)；

13、步骤十：在海域全局坐标系x下对潜艇及围绕潜艇随动流场域网格模型进行整体位移，其中，整体位移矢量dsubmarine＝v•θt；

14、步骤十一：在船体坐标系xsubmarine_local下对潜艇及围绕潜艇随动流场域网格模型进行整体旋转，其中，整体旋转矢量θsubmarine＝ω·δt；其中ω为潜艇角速度矢量；

15、步骤十二：对流场域网格模型基于重叠网格理论进行挖洞处理，并在流场域网格模型和潜艇及围绕潜艇随动流场域网格模型各自重叠网格边界单元处，建立两者间流场信息插值，同时以重叠网格挖洞后潜艇及围绕潜艇随动流场域网格模型的重叠边界单元尺寸为基准，对流场域网格模型重叠边界单元及其周边三层单元进行局部加密；

16、步骤十三：对流场域网格模型中位于海水分界面附近的四层网格重新进行局部加密；其中，海水分界面为淡海水和浓海水的分界面，沿海洋深度方向将海水由上至下根据设定海水浓度阈值的分为淡海水和浓海水两层；

17、步骤十四：根据当前局部坐标系xfluid_local原点所对应的海床信息，对内孤立波波速c、海水分界面到当前海床距离h2和波形方程η(x′,t)进行更新，对步骤三所述速度入口型边界流场速度进行更新，其中，t＝n·δt为当前计算时间，x′为流场域网格模型局部坐标系xfluid_local下总体装配网格模型节点x轴向坐标；

18、步骤十五：对所述总体装配网格模型进行流场计算；

19、步骤十六：计算潜艇受力矢量f＝σepe·ne-∑eτe·ae和潜艇力矩矢量m＝∑e[de×(pe·ne)]-∑e[de×(τe·ae)]，其中，pe为潜艇壁面单元e表面流场压力，ne为潜艇壁面单元e法线方向矢量，τe为潜艇模型壁面单元e所受剪切应力，ae为潜艇模型壁面单元e面积，de为潜艇重心到该壁面单元e距离矢量，∑e表示对潜艇表面所有单元求和；

20、步骤十七：求解方程和更新潜艇质心运动速度矢量v和潜艇角速度矢量ω，其中，m为当前潜艇排水量，i为当前潜艇惯性张量；

21、步骤十八：对步骤十所述整体位移矢量dsubmarine和步骤十一所述旋转矢量θsubmarine进行更新，并输出；

22、步骤十九：令n＝n+1，若n·δt＞tend，则结束计算，否则返回步骤九，其中，tend为计算结束时间。

23、作为本发明的进一步技术方案，步骤一所述流场域几何模型为立方体。

24、作为本发明的进一步技术方案，步骤一所述流场域几何模型位于计算海域的海面和海床之间，与海面平行，且海水分界面贯穿流场域几何模型。

25、作为本发明的进一步技术方案，步骤八所述初始化具体步骤为：

26、步骤1：根据计算海域海床平均深度信息以及海水分界面深度，确定内孤立波理论模型；

27、步骤2：根据内孤立波理论模型确定所述波形方程η(x′,t)，对步骤八中的总体装配网格模型基于vof方法进行流体体积分数初始化，其中，所述波形方程η(x′,t)中的t初始取值为0；

28、步骤3：计算流场域网格模型上表面边界压力pup；

29、步骤4：以流场域网格模型上表面边界压力pup为基准，采用静水压力理论对步骤三所述压力出口型边界流场压力和步骤八所述总体装配网格模型流场压力进行初始化；

30、步骤5：在流场域网格模型局部坐标系xfluid_local下对步骤三所述速度入口型边界流场速度和步骤八所述总体装配网格模型流场速度进行初始化；

31、步骤6：对流场域网格模型基于重叠网格理论进行挖洞处理，并在流场域网格模型和潜艇及围绕潜艇随动流场域网格模型各自重叠网格边界单元处，建立两者间流场信息插值，同时以重叠网格挖洞后潜艇及围绕潜艇随动流场域网格模型的重叠边界单元尺寸为基准，对流场域网格模型重叠边界单元及其周边三层单元进行局部加密；

32、步骤7：对流场域网格模型海水分界面附近四层网格重新进行局部加密。

33、作为本发明的进一步技术方案，步骤3中所述流场域网格模型上表面边界压力pup＝pair+ρ1·g·depup，其中，pair为海面空气压力，ρ1为淡海水密度，g为重力加速度，depup为流场域网格模型上表面边界所处深度。

34、作为本发明的进一步技术方案，步骤十五中对所述总体装配网格模型采用非定常simple算法进行流场计算。

35、作为本发明的进一步技术方案，步骤十五所述非定常simple算法的控制方程中的速度变量采用相对速度vrelative＝vlagrange-vfluid_mesh代替，其中，vlagrange为基于拉氏描述的流体材料运动速度，vfluid_mesh为矢量(c,0,vsubmarine_z)在海域全局坐标系x下的投影矢量。

36、作为本发明的进一步技术方案，步骤十四中，所述速度入口型边界处流体若为淡海水，则其在流场域网格模型局部坐标系xfluid_local下速度矢量为若为浓海水，则其在流场域网格模型局部坐标系xfluid_local下速度矢量为h1为海水分界面到海面距离，h2为海水分界面到海床距离。

37、作为本发明的进一步技术方案，步骤一、步骤四中，对所述流场域几何模型、所述潜艇及围绕潜艇的随动流场域几何模型均采用笛卡尔网络进行网格划分。

38、作为本发明的进一步技术方案，步骤5中，淡海水和浓海水在流场域网格模型局部坐标系xfluid_local下速度矢量分别为和其中t取值为0，h2为海水分界面到海床距离。

39、本发明采用以上技术方案与现有技术相比，具有以下技术效果：

40、本发明采用局部网格加密、运动流场域方法，在大幅度降低网格数量的同时，大幅提高关键流场区域细节描述，充分考虑海床对内孤立波的影响；

41、本发明采用重叠网格技术降低网格变形造成计算终止的风险；

42、本发明公开的技术无需考虑传统造波技术产生的数值误差，可实现大尺度空间域和时间域内孤立波影响评估，与现有计算方法相比，具有计算量小，成功率高，造波简单，数据易采集，经济性高等特点。