一种获取水下航行器直航附加质量方法与流程

1.本发明属于船舶与海洋工程技术领域，具体涉及一种水下航行器直航附件质量的求取方法。


背景技术：

2.水下自主航行器(autonomous underwater vehicle，auv)是一种常用于海洋探测、监察的非标设备。auv会根据用户使用需求而对壳体进行定制，因此auv外形不一。为了更好地操纵auv，就必须先了解其在水中的操纵性。附加质量是物体在理想流体中作非定常运动是所受到的水动力，是潜艇操纵性的其中一个参数。
3.获得非标auv附加质量的方法有efd(empirical fluid dynamics，经验流体力学),cfd(computational fluid dynamics，计算流体动力学)、hess
‑
smith面元法和水下实验。其中水下实验是最准确的方法，但成本高、耗时长。而efd只能计算标准形状的附加质量。hess
‑
smith面元法计算附加质量的精度满足工程应用需求，但程序代码复杂，较难快速上手。
4.ansys是成熟的流体cfd仿真软件，能用于计算附加质量，常见的利用ansys fluent计算附加质量使用了动网格技术，分别模拟物体在水中做规律的正弦运动和水平匀加速直线运动，并通过分析仿真所得数据求得直航方向附加质量。由于需要使用动网格技术，网格质量对仿真十分关键。使用动网格仿真，势必会造成部分网格移动重建，且为保证cfd仿真精度，往往越复杂外形的物体表面网格尺寸越小，小尺寸网格限制了动网格仿真的单位时间步长，从而降低了求解速度而无法快速求解附加质量。


技术实现要素：

5.本发明提出一种基于ansys的cfd求解外形不一的auv直航方向的附加质量的方法，与传统的基于动网格的cfd方法不同，该方法不需要使用动网格，在保证计算精度的同时降低了网格要求且极大的提升了计算速度，其采用如下技术方案予以实现：
6.一种获取水下航行器直航附加质量方法，包括：
7.a、auv模型导入；
8.b、网格划分；
9.c、网格独立性验证，获得较优网络；
10.d、进行时间步长独立性验证，获得较优的瞬态仿真时间步长；
11.e、网格导入；
12.f、将udf程序编译并导入fluent中进行瞬态仿真
13.auv做匀速直线运动时，auv的推进力f
p
大小等于匀速运动时的阻力drag；而当auv做加速恒为a的水平直线运动时，将匀加速运动过程离散化，选取vt前后时间段的推力值，使用最小二乘法拟合t时刻前后一段时间的数值，来推算某一时刻的推力f
vt
[0014][0015]
结合公式f
vt
‑
drag
vt
‑
ma＝m
11
·
a得出：
[0016]
其中ε为修正系数。
[0017]
进一步地，所述步骤b包括：创建外流域和内流域，内流域划分成非结构网格，外流域为划分成结构网格，同时在auv表面添加边界层网格。
[0018]
进一步地，基于步骤b中的网格，生成一系列连续细化的网格，网格数量变化比率为
[0019]
与现有技术相比，本发明的优点和积极效果在于：
[0020]
本发明与传统的基于动网格的cfd方法不同，该方法不需要使用动网格，通过改变流体速度模拟auv变速的过程，监控阻力值并通过计算求解附加质量。本发明在保证计算精度的同时降低了网格要求且极大的提升了计算速度。为验证本发明准确性，对标准椭圆球分别进行cfd仿真及理论计算。本发明得到的附加质量值与理论计算值偏差小于3％，说明在设计初期，用该方法计算任意外形的auv直航方向的附加质量结果可靠。
附图说明
[0021]
图1为本发明实施例获取水下航行器直航附加质量方法流程图；
[0022]
图2为本发明实施例auv坐标系；
[0023]
图3为本发明实施例椭圆球模型；
[0024]
图4为本发明实施例流域设置；；
[0025]
图5为本发明实施例网格整体图及截面图；
[0026]
图6为本发明实施例suboff模型；
[0027]
图7为本发明实施例阻力与时间关系图一；
[0028]
图8为本发明实施例阻力与时间关系图二；
[0029]
图9为本发明实施例阻力与时间关系图三。
具体实施方式
[0030]
ansys fluent是成熟的流体仿真软件，能准确地对不同形状的物体在不同的工况下求解rans方程。mochammad、raza、姚保太等人利用动网格，分别模拟椭圆球在水中做规律的正弦运动和水平匀加速直线运动，并通过分析仿真所得数据，准确地求得直航方向附加质量。由于需要使用动网格技术，网格质量对仿真十分关键。仿真过程中部分网格会移动重建，所以动网格仿真的最大步长取决于最小网格尺寸，对于复杂外形的物体便无法快速求解附加质量。本实施例提出一种无需使用动网格快速求解直航方向附加质量的cfd求解方法。
[0031]
为了能够更加清楚地理解本发明的上述目的、特征和优点，下面结合附图及实施例对本发明做进一步说明。
[0032]
参考图1，本实施例基于ansys fluent的cfd求解外形不一的auv直航方向的附加质量的过程如下：
[0033]
步骤a、导入auv三维模型于ansys icem cfd软件中
[0034]
当然也可选择将三维模型导入至其他网格划分软件中，如gambit、meshing等。
[0035]
步骤b、网格划分
[0036]
实现过程：创建外流域和内流域，内流域划分成非结构网格(非结构网格优点：适应性极强，能较好地处理复杂平面。)外流域为划分成结构网格(结构网格的优点在于：占用计算资源少；计算速度快)。同时在水下航行器auv表面添加边界层网格(边界层是高雷诺数绕流中紧贴物面的粘性力不可忽略的流动薄层。流体流经固体壁面时会在紧贴壁面的区域内形成边界层，边界层内流体具有很大的速度梯度。因此本实施例对水下机器人进行流体分析时考虑边界层，添加边界层网格可以在保证计算精度的同时极大地降低网格数量以提升仿真运算速度和减低计算机资源使用。)
[0037]
作为替代方案，也可只创建一个流域，将流域内网格全部划分为结构网格或全部划分为非结构网格。
[0038]
步骤c、进行网格独立性验证，得到较优网格
[0039]
一般而言，网格数量越大，运算精度越准，但运算时间越久，本实施例选取的较优网格就是网格数量尽量小的同时，运算精度也比较好。为合理地使用计算机资源，较快地完成瞬态仿真，本实施例基于步骤b初始网格，生成一系列连续细化的网格，网格数量变化比率为通过流体仿真选择网格数量较少且仿真数据变化不大的网格。如表2所示，mesh2,mesh 3,mesh 4仿真计算结果相似，所以选取数量最少的mesh2为较优网格。
[0040]
步骤d、进行时间步长独立性验证，得到较优的瞬态仿真时间步长
[0041]
与较优网格相似，本实施例选取的较优时间步长就是时间步长尽量大的同时，运算精度也比较好。
[0042]
步骤e、选定一个直航速度
[0043]
将较优网格导入ansys fluent中，并对其进行直航方向1m/s匀速直线运动的流体仿真，并储存数据drag
1m/s
，水下航行器航速设为0
‑
2.5m/s本实施例优选1m/s。
[0044]
步骤f、将udf程序(本实施例优选加速度a为0.1m/s2为示例)编译并导入fluent中进行瞬态仿真
[0045]
udf规定：0到2.4s入口速度为0.8m/s，2.4s后入口速度为[0.8+a
·
(t
‑
2.4)]m/s，运算时间t设置为0
‑
5s，并记录每个时间步长的阻力(推力)数据。本实施例优选加速度a为0.1m/s2，可替换加速度a数值大小，但加速度不能过高；初始速度0.8m/s可替换，速度大小尽量接近第五步的直航速度；初始速度仿真时间2.4s可替换，考虑到瞬态仿真初期仿真会产生较大的偏差，先使用恒定来流速度仿真至匀速直线运动阻力基本保持不变后再进入匀加速直线运动，要求仿真时间后在瞬态匀速仿真基本稳定。
[0046]
在水下直航时，auv只受到由加速度引起的惯性力和速度引起的粘性力。当auv做匀速直线运动时，auv的推进力f
p
大小如下所示，其中drag为匀速运动时的阻力。
[0047]
f
p
＝drag
[0048]
而当auv做加速恒为a的水平直线运动时，在t时刻auv的瞬时推力大小等于在t时刻由速度产生的阻力大小加上在匀加速度a下的auv质量与附加质量所产生的惯性力，关系
式如下：
[0049]
f
pt
‑
drag
t
‑
ma＝m
11
·
a
[0050]
将匀加速运动过程离散化，则匀加速运动近似为auv在时间差δt内做匀速直线运动，并将上一迭代的结果带入下一迭代的开始。利用cfd瞬态计算，通过监控auv的阻力实现记录每个时间差的阻力drag
t
‑
nδt+i
。并通过最小二乘法拟合t时刻前后一段时间的数值，来推算某一时刻的推力f
pt
，以保证拟合而得的推力f
pt
为t时刻的瞬时推力。
[0051][0052]
将如上三式结合，便可得：
[0053]
其中ε为修正系数。
[0054]
为验证本发明提出的方法能准确计算水下航行器(auv坐标系参考图2)的附加质量，选取3个简单的、可进行理论计算的椭圆球模型，如图3所示，l为椭圆球长度，d为椭圆球截面直径，m为椭圆球的质量。
[0055]
为降低尾流效应和增加水动力性能，auv的长径比一般为10。实际生产的auv长径比一般为7
‑
10。本实施例选取长径比为5、7、10的3个椭圆模型进行仿真验证。模型的具体几何参数如表1所示。
[0056]
表1椭圆球模型参数
[0057][0058]
一个任意形状的物体运动时共有36个附加质量，可列成方阵如下：
[0059][0060]
由势流理论式(i,j＝1，2，
…
，6)，可得m
ij
＝m
ji
，式中s为物体表面积；n为微元面积ds的外法线方向；分别是物体延坐标轴x、y、z以单位速度平移运动所引起的流体速度势；分别是物体以单位角速度绕x、y、z轴作转动运动所引起的流体速度势。由于auv基于xoz面对称，且各附体较小，对附加质量影响较小，可近似看作auv基于xoy对称，所以方阵可简化为如下：
[0061][0062]
方阵中除了m
11
和m
44
，其他均可通过平面运动机构试验模拟(planar motion mechanism，pmm)纯升沉、纯横荡、纯俯仰、纯艏摇直接或间接求得。
[0063]
blevins给出了椭圆球的附加质量m
11
的经验公式：
[0064][0065]
式中ρ是水的密度，α是由blevins测量的经验参数，由椭圆的长度与直径之比确定，α值如表1所示，l为长轴长度，d为短轴长度。由公式(5)计算得椭圆1，椭圆2和椭圆3的m
11
理论值如表1所示。
[0066]
仿真模型的流域尺寸如图4所示：速度入口边界距离椭圆前端13倍椭圆最大截面直径，自由出口距离椭圆尾端33倍椭圆最大截面直径。外流域采用占用计算资源少的结构网格，内流域采用非结构网格。
[0067]
如图5所示，边界条件设置如下：入口设置为速度入口(velocity
‑
inlet)，出口设置为自由流(outflow)，wall设置为对称面(symmetry)，将流域交界面interface将结构网格与非结构网格的节点融合并将类型设置为interior。
[0068]
对不可压湍流模型仿真，使用fluent进行有限体积法求解雷诺平均方程，采用有限体积法对二阶迎风格式的控制方程进行离散。压力
‑
速度方程采用半隐式方程组(simple)，湍流模型选择自由剪切流动和壁面限制边界层流动优势明显的sst k
‑
omega，同时需要在网格模型中添加边界层。边界层的首层高度为:
[0069][0070]
其中y
+
是无量纲，υ为流体运动粘度；u
*
为近壁面摩擦系数，τ
w
为壁面剪切应力，本实施例取y
+
＝30。
[0071]
在cfd分析之前，需网格独立性试验，选取合适的网格文件进行后续分析。以model3为对象，来流速度为1m/s，网格独立性试验如表2所示。
[0072]
表2网格独立性验证
[0073]
[0074]
本实施例选取mesh2作为后续model3的仿真对象，同时按照相同参数绘制model1和model2的网格文件用于后续仿真。
[0075]
为验证本实施例直航仿真的准确性，本实施例对ittc公布的suboff标准模型进行不同航速下的直航仿真，suboff标准模型主要几何参数如表3所示；
[0076]
表3 suboff物理参数
[0077][0078]
由于suboff模型与椭圆3的外形尺寸相似，按照mesh2相同网格参数对suboff标准模型进行三维结构网格划分，suboff模型如图6所示：
[0079]
通过cfd计算得出的阻力值与ittc实验值进行对比，如表4所示，仿真值与实验值偏差不超过5.0％，可以验证本实施例的仿真方法可靠。
[0080]
表4仿真准确性验证
[0081][0082][0083]
将model1，2，3的网格分别导入fluent中求解，选用稳态求解器，求解3个模型分别在水下1m/s直航时所受阻力drag
1m/s
。记录后选用瞬态求解器，导入udf，udf规定：0到2.4s入口速度为0.8m/s，2.4s后入口速度为[0.8+a
·
(t
‑
2.4)]m/s。运算2.4s后停止运算。以model1为模型，加速度为0.1m/s2为例，监控瞬态每个时间步长的阻力，阻力与时间关系图如图7
‑
图9所示。
[0084]
代入式求解。同理将model1，model2，model3分别以加速度0.1m/s2，0.2m/s2，0.4m/s2进行数值计算及数值拟合，结果如表5所示。
[0085]
表5仿真结果
[0086][0087]
本发明得到的附加质量值与理论计算值偏差小于3％，说明在设计初期，用该方法计算任意外形的auv直航方向的附加质量结果可靠。
[0088]
以上所述，仅是本发明的较佳实施例而已，并非是对本发明作其它形式的限制，任何熟悉本专业的技术人员可能利用上述揭示的技术内容加以变更或改型为等同变化的等效实施例应用于其它领域，但是凡是未脱离本发明技术方案内容，依据本发明的技术实质对以上实施例所作的任何简单修改、等同变化与改型，仍属于本发明技术方案的保护范围。