涡流场中波动鱼形体推进特性数值模拟

李永成，陈纪军，闫 朋

（中国船舶科学研究中心，江苏 无锡 214082）

0 引 言

自然界中的鱼类经过长时间的进化之后已具备非凡的运动能力，其运动方式主要分为鳗鲡模式、鲹蝌模式和月牙尾推进模式等3种。在鲹蝌模式和月牙尾推进模式中，鱼类仅部分身体参与摆动，其余部位保持刚性，典型代表有鲸鱼等。这2种运动方式可抽象为拍动运动，通常用耦合的升沉运动和俯仰振荡运动代替。在鳗鲡模式中，鱼类全身作横向摆动，以行波方式向前运动，典型代表有鳗鱼等。这种运动方式可抽象为鱼形体波动运动，在理论分析或数值计算中，通常用平板或对称翼型代替鱼身，沿身长方向的各点作横向振动。对于横向振动的具体运动形式，一种较为常用的描述方式是平板或对称翼沿身长方向作行波运动，身长方向上各点的横向位移可通过正弦函数或余弦函数表示。

国内外已有很多学者对鱼形体波动推进进行研究。

1) 在理论研究方面：LIGHTHILL提出了细长体理论，用以求解波动鱼形体的受力问题；WU提出了三维波动板理论，成功地解决了无黏流场中三维波动鱼形体的数值计算问题；CHENG等提出了半解析半数值的涡格法，该方法已成功应用于三维波动鱼形体流场模拟中。

2) 在试验研究方面：BARRETT等通过研制仿生波动机器鱼，证实了与相同尺寸的刚性鱼相比，波动机器鱼达到相同航行姿态耗费的能量更少；LAUDER应用DPIV（Digital Particle Image Velocimetry）研究了波动鱼形体长鳍的漩涡尾迹与尾鳍的漩涡尾迹的相互影响，结果表明，长鳍与尾鳍相互作用的水动力学效应可使鱼形体的推进力得到增强；KATO对低速条件下鱼形体波动运动特性进行了分析，其研制的第四代样机具有良好的操作性。

3) 在数值计算方面：EPSTEIN等采用RANS方法将鱼形体简化为正弦行波运动，获得了其推力与行波运动参数之间的关系；RAHMAN等以乌贼为仿生对象，将波动鳍简化为等摆幅波动的矩形鳍面，采用基于泊松方程控制生成的体贴动网格，计算得到了鳍面上的压力分布特征；潘定一等和LI等采用浸没边界法对二维波动、三维波动和形体进行了低雷诺数下的直接数值模拟研究，确立了最佳运动参数。

在上述关于鱼形体波动的研究中，一般将鱼形体置于均匀来流或静止流体中，流场可认为是稳态流场。然而，在真实的自然界中存在着各种障碍物（如珊瑚和礁石等），流场中存在着各种漩涡鱼形体，其周围的流场通常较为复杂。鉴于此，本文主要研究涡流场（非稳态流场）中波动鱼形体的推进特性。通过在鱼形体上布置一个静止方柱模拟非稳态流场，系统分析静止方柱与鱼形体之间的距离和斯特哈尔数对波动鱼形体推进特性的影响，为以后小型仿生潜器研制提供参考。

1 计算模型

图1为物理模型示意。计算模型为二维NACA0012剖面翼型，弦长为

C

，与一个方柱前后串联布置，翼型顶端与方柱后端之间的距离为

S

，方柱的边长为

a

。选取一个矩形计算域，来流从左至右依次流过方柱和鱼形体，以翼型弦长为特征长度，整个计算域的范围在流向上为26倍弦长，在垂直方向上为16倍弦长。关于边界条件设置，除了将计算域右侧设为压力出口边界条件，计算域其他边界均设为速度入口边界条件。

图1 物理模型示意

鱼形体的头部被固定，身体的其他部分作行波运动，变形起始位置与鱼形体头部之间的距离为0.3

C

。对于作行波运动的鱼形体而言，其横向位移的表达式为

式(1)和式(2)中：

A

(

x

)为幅值函数；

a

、

a

和

a

为用于控制鱼形体的波动幅值，根据参考文献的建议，

a

一般取为－0.02；

k

为波数，其表达式为

k

＝2π/

λ

，

λ

为波长；

f

为运动频率；

x

≥0.3

C

。根据以往的经验，为获得最佳的推进效果，

a

与

a

的关系为

a

＝－2

a

。图2为网格划分示意，考虑到翼型一直都在作横向波动，除了翼型周围采用非结构化网格，计算域其他地方均采用结构化网格。翼型表面第一层网格的高度满足

y

≤1，使得第一层网格位于黏性底层区。

图2 网格划分示意

2 数值方法及验证

2.1 控制方程

考虑不可压缩流动湍流问题，其运动控制方程为

式(3)和式(4)中：

u

(

i

＝ 1,2)为流体运动速度；

x

(

i

＝ 1,2)为空间坐标；

P

为流体压强；

γ

为运动黏性系数；

γ

＝

ck

ε

为湍流黏性系数，其中，

k

为湍流动能，

ε

为湍动能耗散率，

c

为常数。

2.2 动网格技术及精度验证

在鱼形体波动推进过程中，翼型边界始终处于运动状态，属于典型的动边界问题。对于动边界问题的处理，采用商用软件Fluent中的动网格技术可有效捕捉运动的边界，通过用户自定义函数（User Defined Function, UDF）控制边界运动即可实现预定的运动。因此，本文采用Fluent中的动网格技术进行数值求解。

Fluent动网格技术中主要有3种调节计算内部网格节点的算法，分别为弹性光顺、局部重构和铺层。在弹性光顺模型中，网格的边被理想为网格节点间相互连接的弹簧，通过胡克定律计算网格边界运动对整个弹簧系统受力的影响，进而通过调整各网格节点的位置，使之达到新的平衡状态。该方法的网格变形较小，计算效率较高，尤其适于求解低频率的变形运动问题。因此，本文选择弹性光顺模型处理运动边界问题。

为验证本文计算方法的可行性，对一个作行波推进运动的水翼进行数值模拟，计算模型采用与参考文献[16]相同的模型，即弦长

c

＝0.1m的二维NACA0012水翼，来流速度设为

v

＝ 1.2m/s ，其他运动参数设置：

a

＝ 0.02；

a

＝－0 .08；

a

＝ 0.16；

λ

＝1.3；

f

的变化范围为0.6～2.1Hz。计算摆动水翼在不同运动频率下的推力系数

C

和推进效率

η

，计算结果见图3，具体表达式为

图3 数值方法验证结果

式(5)中：

T

为变形翼受到的推力；

P

为系统输出功率；

P

为系统输入功率，其含义为鱼形体横向振动时克服流体阻力所做的功；

p

和

v

＇分别为变形翼表面压力和速度分布；

Γ

为变形翼边界；

t

为某一运动时刻；

ρ

为流体密度；

C

为水翼的弦长；

l

为水翼的周长。

由图3可知，本文的计算结果与参考文献[16]中的计算结果能较好地吻合，因此本文建立的数值计算方法可用来研究摆动可变形翼的仿生推进性能。

3 计算结果和讨论

3.1 涡流场与均匀流场中波动鱼形体推进特性比较

首先对涡流场中波动鱼形体与均匀流场中波动鱼形体的推进特性进行比较。计算工况见表1，其中：SF1～SF4为均匀流场工况；UF1～UF4为涡流场工况。不同工况组合下鱼形体推力系数和系统输入功率随运动频率的变化曲线见图4。

表1 计算工况

图4 不同工况下鱼形体推力系数和系统输入功率随运动频率的变化曲线

1) 由图4a可知，对于推力系数，在均匀流场中作波动推进的鱼形体的推力系数小于在涡流场中作波动推进的鱼形体的推力系数，尤其是在高运动频率下，二者的差异更为显著。此外，在低频率运动工况下，均匀流中的波动鱼形体获得的推力为负值，这意味着鱼形体无法产生向前的推力。相比之下，涡流场中的鱼形体获得的推力为正值，这意味着涡流场中的卡门涡街对鱼形体的推进有积极作用。

2) 由图4b可知，对于输入功率，相同运动频率下涡流场中的波动鱼形体所需的系统输入功率远低于均匀流场中的波动鱼形体所需的系统输入功率，且运动频率越高，二者的差距越明显。特别是当运动频率

f

=1.5Hz时，涡流场中鱼形体的系统输入功率仅为均匀流场中的1/8。这进一步说明，涡流场中存在的漩涡对鱼形体推进有显著的促进作用。图5为SF和UF工况下鱼形体周围流线和压力（

C

）云图（

t

= 5

T

）。由图5可知：流场上游方柱产生的漩涡流经过其后方的鱼形体时，鱼形体表面的流线发生改变甚至反转，导致克服流体阻力消耗的能量减少；与均匀流场相比，涡流场中鱼形体头部的压力相对较小，这正是涡流场中的鱼形体能获得较强推力的原因。

图5 SF和UF工况下鱼形体周围流线和压力云图

3.2 方柱与鱼形体的间距对鱼形体推进特性的影响

系统分析方柱与鱼形体的间距

S

对鱼形体推进特性的影响。分别选取

S

的值为0.15、0.25、0.50、0.80和1.20，探究不同

S

下鱼形体推力系数和推进效率的变化情况，计算结果见图6。为更方便地对鱼形体的推力系数进行分析，将推力系数分解为压差造成的推力系数

C

和黏性导致的推力系数 v

TC

，即

式(6)中：

T

和

T

分别为压差阻力和黏性阻力。为更直观地显示方柱与鱼形体的相互作用，绘制不同

S

下鱼形体后方的涡流场结构图，结果见图7，其中

Z

为流场

z

方向的涡量，其表达式为

式(7)中：

v

和

v

分别为鱼形体在水平方向和垂直方向上的运动速度。

由图6a和图7可知，涡流场中波动鱼形体受力可分为3种情况。

图7 不同间距下鱼形体后方涡流场结构

1) 当

S

＜0.5

C

时，由于鱼形体与方柱之间的距离很短，鱼形体的存在抑制了方柱的涡环脱落，且鱼形体表面被方向相反的涡环覆盖。鱼形体的波动运动抑制了静止方柱尾流场中剪切层的卷起和漩涡脱落，延伸了吸力区，使得其推力特性增强。此时，黏性阻力产生正的推力系数，且随着

S

的增大，对应的

C

逐渐减小，使得整体的推力系数略有减小。2) 当0.5

C

≤

S

≤1.0

C

时，随着鱼形体与方柱之间的距离逐渐增大，方柱后方逐渐形成完整的上下交错的涡环。鱼形体的推力系数随着

S

的增大而逐渐减小，且黏性阻力产生的推力 v

TC

也逐渐减小，直至为0。3) 当

S

＞1.0

C

时，鱼形体与方柱之间的距离很大，方柱后方脱落的涡环在到达鱼形体头部过程中逐渐耗散，方柱与鱼形体之间距离的影响明显减弱。鱼形体的推力系数迅速减小，黏性阻力产生的推力 v

TC

为负值。由图6b可知，涡流场中鱼形体的推进效率随着

S

的增大呈逐渐减小的趋势，尤其是当

S

＞0.5

C

之后，鱼形体的推进效率更是大幅下降，这与推力系数的变化是一致的；当

S

继续增大时，鱼形体与方柱的相互作用逐渐减弱，推进效率值近乎不再发生变化。

图6 平均推进力系数和推进效率随鱼形体与方柱间距的变化曲线

4 结 语

本文通过在鱼形体上游布置一个静止方柱模拟涡流场环境，采用动网格技术对涡流场中波动鱼形体的推进特性进行数值模拟研究，系统分析方柱与鱼形体间距对鱼形体推进特性的影响。数值计算结果表明：与均匀流场中的波动鱼形体相比，涡流场中的波动鱼形体可获得更大的推力，且消耗的功率较少；涡流场中波动鱼形体的受力情况随着方柱与鱼形体间距的变化而有所不同，随着二者间距的增大，鱼形体的推力呈逐渐减小的趋势。