基于柔性可折叠翼的水下滑翔机建模与运动特性分析

李振宇，陈质二，俞建成，杨 庆，4

（1.中国科学院沈阳自动化研究所机器人学国家重点实验室，沈阳 110016；2.中国科学院机器人与智能制造创新研究院，沈阳 110169；3.中国科学院大学，北京 100049；4.东北大学机械工程与自动化学院，沈阳 110819）

0 引 言

水下滑翔机由于自身航速较低且续航时间有限[1]，其通常由母船送至目标海域布放。在多数情况下观测海区距离陆基较远，航渡时间较长（一般需要数天），尤其是在面对海上突发事件和海洋现象时滑翔机的响应速度较慢，同时对于具有隐蔽布放需求的场合，如图1所示，由于母船形体过大，布放活动易于暴露。因此通过载机空投和潜艇鱼雷管发射来实现水下滑翔机的部署逐渐得到重视和研究。相比于传统的刚性翼面水下滑翔机，柔性可折叠翼水下滑翔机的机翼在入水前折叠，入水后自动弹开，可以解决传统滑翔机空投时刚性翼板存在的入水冲击损坏[2]和潜艇水下发射管径受限问题。

近年来，柔性可折叠翼较多地应用在航空信息以及仿生领域，尤其是扑翼类微型飞行器方面。韩国建国大学的Muhammad 等[3]根据英国埃克塞特大学的Haas 等提出的四平面理论，利用形状记忆合金丝作为驱动器，成功实现了可折叠翼的自动折叠展开动作；后来Truong 等[4]提出了基于四杆机构的可折叠翼，该翼最大优势在于结构简单，缺点在于需要依靠外力推动连杆结构才能实现折叠，翼在扑动中形态发生变化后无法及时修正；南京航空航天大学顾鑫[5]对柔性折叠翼飞行器飞行动力学进行研究，分析了弹性变形和折叠角速度对纵向运动参数的影响规律；斯坦福大学的Stowers等[6]仿照鸟类和蝙蝠的翅膀形态，设计了一种利用离心加速度使机翼从初始折叠状态转变为展开状态的仿生可折叠翼；浙江大学李铁风等[7]研制出仿生软体柔性机器鱼，其由软体人工肌肉驱动一对翼状柔性胸鳍，通过节律性扑翅实现游动，翼翅选用耐高压的软体材料，可以承受万米级别深海静水压力。

国内外学者也对具有变翼功能的水下滑翔机开展了研究。大阪大学Arima 等[8-9]设计出可转动翼板翻转角度的“ALEX”滑翔机，其通过建立滑翔机的动力学模型并对不同翼板翻转角度的滑翔机进行运动仿真，仿真结果表明，可转动翼板角度的滑翔机具有更小的滑翔攻角和更好的运动性能；天津大学杨志金等[10]提出一种具有变体功能的可变翼混合驱动水下滑翔机，通过改变机翼展弦比、后掠角和翻转角可提高滑翔机在不同工况下的航行性能；西北工业大学田文龙等[11]对可变翼混合驱动水下滑翔机进行研究，基于Newton-Euler方程建立了六自由度运动数学模型，通过求解稳态滑翔运动平衡方程，定量分析了可控翼翻转角对滑翔机滑翔性能的影响。上述多为机翼自由度数目增加的变翼研究，水下滑翔机通过载机空投和潜艇鱼雷管发射时需要机翼能够完全折叠。

本文以中国科学院沈阳自动化研究所自主研制的海翼1000mini水下滑翔机为柔性可折叠翼搭载基体，研究水下滑翔翼面展开后的水动力特性并建立其在流体中的动力学模型。根据水动力与运动仿真实验的结果，分析翼形变与水下滑翔机水动力特性和运动参数之间的关系，可为后续滑翔机控制算法设计和运动优化提供理论依据。

1 水动力参数辨识

1.1 翼面形位描述与形变量定义

柔性可折叠翼水下滑翔机机翼展开后柔性翼面受周围流场的扰动会产生一定的形变，折叠翼的翼面支撑杆为不锈钢材质，相较于翼面可视其为刚性杆。根据刚柔耦合多体系统节点描述方法[12]，建立如图2 所示的惯性坐标系（I-i,j,k）和机翼坐标系（S-x,y,z），翼面上某一点f0的形变为机翼坐标系相对于惯性坐标系的位置矢量和柔性翼面相对于机翼坐标系变形矢量的叠加，对于翼面上某一点fi变形后对应的Fi在惯性坐标系下表示为

式中，m为Fi在惯性系中位置向量，K为机翼坐标系S在惯性坐标系I中对应的位置向量，A为方向余弦矩阵，li为翼未发生形变时fi点在机翼坐标系中的位置向量，hi为柔性翼的相对变形量。

在惯性系中的机翼坐标系的位置向量M与方向余弦矩阵A均为定值，所以将形变描述简化至机翼坐标系中，下文中均选定翼面外边缘弧顶至翼面坐标系y轴的标量距离h=‖h‖作为翼面的形变量。

1.2 水动力修正计算公式引入

水下滑翔机航行时所受到的水动力包括滑翔机的升力L、阻力D与俯仰力矩M，当滑翔机以小攻角航行，不考虑翼形变时滑翔机的升力、阻力和力矩满足如下关系式[13]：

式中，ρ为海水密度，A为横截面积，V为航行速度，α为航行攻角，CD为阻力系数，CL为升力系数，CM为俯仰力矩系数。KD0、KL0、KM0为与滑翔机直航特性有关的水动力系数，KD、KL、KM为与滑翔机攻角特性有关的水动力系数。

翼形变h会对滑翔机升力、阻力和俯仰力矩产生影响，可将式（2）中水动力系数写为与翼形变有关的变量，因此将柔性可折叠翼水下滑翔机水动力计算公式修正如下：

式 中，KWD0(h)、KWD(h)、KWL0(h)、KWL(h)、KWM0(h)、KWM(h)为修正后的与滑翔机翼面形变有关的水动力系数。

1.3 柔性可折叠翼结构设计与工作原理

如图3 所示，柔性可折叠水下滑翔机在其艉部密封钢圈和机身回转体上安装一对可折叠机翼。其中，翼面支撑杆底端转动块和基体铰接在一起，铰接轴外围并排安装三个扭簧，扭簧驱动转动块绕铰接轴旋转。支撑杆底端可转动体和基体之间设计一对棘轮自锁机构，使开合后的翼处于目标位置且不会受外界影响产生复位和抖动。滑翔机入水前用水溶性胶带捆绑使其处于收合状态的机翼，滑翔机入水后水溶性胶带遇水溶解断裂，折叠翼迅速自动展开，水溶胶带断裂瞬间翼面展开过程如图4所示。

1.4 水动力仿真实验设计

本文通过STAR-CCM+水动力仿真软件构建仿真模型，计算流域是长为10l，宽和高均为3.5l的长方体（l为滑翔机的长度），如图5所示。滑翔机主轴线与流域的长轴重合，流域速度进口端为长方体的前端面，与滑翔机艏部距离为4l，流域压力出口端为长方体的后端面，与滑翔机艉部距离为5l。计算水域采用多面体网格进行网格划分，其CFD 分析模型的网格离散方式如图6 所示，为提高计算精度，采用棱柱层网格对水下滑翔机近壁面湍流进行捕捉，其中棱柱层网格近壁厚度设置为1 mm，层数设为8 层，网格增长率为1.3。仿真模型中加入面控制对滑翔机的机体和翼的表面进行了加密。建立四组不同网格尺寸的计算模型开展网格无关性验证，四组网格尺寸下网格数量分别为80万、120万、150万、200万，仿真计算均在1.0 m/s航速和4°攻角下进行，不同网格数量下阻力的仿真计算结果如表1所示。通过对比分析可知，当网格数达到150万时滑翔机的阻力基本上不再发生变化，综合考虑计算精度和效率，网格数量选取为150万。

表1 不同网格数量阻力计算结果Tab.1 Computational results of lift-drag ratio of different mesh numbers

已知流体运动粘性系数υ= 1.003 74×10-6/(m2·s-1)，设定滑翔机最大航速为1 m/s，机身总长度为1.08 m，求得对应的雷诺系数为

所以计算过程应选择湍流模型，由于Re在106数量级，本文在Star CCM+仿真软件中选择低雷诺数k-ε湍流模型进行流体仿真计算。

1.5 水动力计算公式的确定

为确定不同翼形变量时柔性可折叠翼水下滑翔机水动力特性变化规律，将不同形变后的翼面简化为刚性翼处理，对0～6 mm 形变时滑翔机的水动力特性开展仿真实验，以0.5 m/s 航速下水动力实验结果为例，水动力数值计算结果如图7所示。

由图可知，翼面存在形变时滑翔机升力、阻力、俯仰力矩与攻角分别满足线性、二次型、线性关系，符合海翼系列水下滑翔机不同攻角下的水动力变化规律[14]。其中，升力和阻力同时随翼面形变的增加而减小，俯仰力矩随翼形变量增加先变大后变小。翼面形变的增加使0°～8°攻角下对应的升阻比逐渐减小，当攻角大于10°时，不同翼面形变下的升阻比基本保持不变。根据图7 中得到水动力实验仿真结果分析式（3）中升力、阻力、俯仰力矩对应的水动力系数KWL0、KWL、KWD0、KWD、KWM0、KWM与形变量之间的关系，各水动力系数拟合计算结果如图8 所示，可以看出KWM0、KWM与形变量近似呈正弦关系，KWD0、KWD、KWL、KWL0与形变量近似呈线性关系，则柔性可折叠翼水下滑翔机水动力计算满足下式：

用( · )表示D、L、M，则( ·)a0-d0、( ·)a-d分别表示水动力计算修正公式中直航特性与攻角特性水动力系数，将拟合后的结果代入上式，得到修正后的柔性可折叠翼滑翔机水动力计算公式为

2 动力学建模

2.1 坐标系的建立与坐标变换

如图9所示，建立惯性坐标系E0-ijk，向量b=[x,y,z]T和向量S=[φ,θ,ψ]T表示滑翔机的浮心在惯性坐标系中的位置和姿态。再以滑翔机浮心为坐标原点建立机体坐标系e0-e1e2e3，分别用V=[u,v,w]T和Ω=[p,q,r]T表示滑翔机在机体坐标系下的速度与角速度。

速度坐标系π0-π1π2π3原点与机体坐标系原点重合，令机体坐标系先绕e2轴转动α，使得e3与π3重合，然后再将旋转之后的坐标系绕π3转动β，此时e1与e2分别与速度坐标系坐标轴π1和π2重合。所以从速度坐标系到机体坐标系之间的转换矩阵可表示为

根据惯性坐标系和机体坐标系之间的转换关系[15]可知

2.2 动力学模型

水下滑翔机是由壳体、可移动质量块、姿态调节系统等组成的多刚体系统[14]。对于图10所示的柔性可折叠翼滑翔机，将其整体分解为三个质点：壳体质点（质量ms，质心rs）、可移动电池块质点（质量mp，质心rp）、浮力调节质点（质量mb，质心rb）。由上述划分可知系统的总质量为

已知电池块质心偏心距为Rp，当其转动角度为γ，沿e1轴移动距离为rp时，其质心rp可表示为

令p和P分别表示滑翔机在惯性坐标系下和机体坐标系下的动量，π和Π表示惯性坐标系和机体坐标系下的角动量，md表示排水量。参照牛顿-欧拉方程，滑翔机系统的动量和受力满足如下关系[13]：

式中，fext和τext分别为外界环境对滑翔机产生的力和力矩，lsE、lpE、lbE分别为图10 中惯性坐标系下壳体质点、电池块质点、浮力质点对应的位置矢量，g表示重力加速度。

定义运算符·为向量“·”的叉乘，根据惯性坐标系和体坐标之间坐标变换关系进一步推导可得柔性可折叠翼滑翔机广义力和力矩在体系中的表达式为

其中，FW=Rfext，TW=Rτext，二者为柔性可折叠翼滑翔机在流体中受到的粘性水动力项，根据式（5）可知，在体坐标系表达式为

式中，V=‖V‖，表示机体坐标系中滑翔机相对于海流的速度，KMR、Kp、Kβ、Kq、KMY、Kr为滑翔机自身的水动力系数。

定义M为系统的广义质量矩阵，其由附加质量矩阵Mf=diag([Mf1Mf2Mf3])、附加转动惯量矩阵IA=diag([If1If2If3])和耦合项组成。将柔性可折叠翼滑翔机机体坐标系下的运动特征用矢量u=[V,Ω]表述，动量特征用η=[P Π]表示，则η和u之间的关系可表示为

上式对时间进行求导，可得

根据上式等可推出u̇的表达式为

再将上式进一步展开得到柔性可折叠翼水下滑翔机运动方程为

3 仿真与实验分析

根据滑翔机动力学模型编写仿真程序并使用四阶龙格库塔法对式（20）微分方程组进行求解[17-18]，得到滑翔机滑翔运动时各运动参量随时间的变化。仿真分析时滑翔机的水动力系数和设计参数如表2～3所示。

表2 柔性可折叠水下滑翔机水动力系数Tab.2 Hydrodynamic coefficients of the flexible foldable-wings glider

表3 柔性可折叠翼水下滑翔机设计参数Tab.3 Design parameters of the flexible foldable-wings glider

3.1 仿真与海试比较分析

为验证本文计算和仿真的正确性，在输入相同控制量条件下，将柔性可折叠翼水下滑翔机0形变仿真结果与海翼1000mini水下滑翔机2020年12月南海海试实验数据进行比对分析，滑翔机的深度和俯仰角的变化如图11 所示，海试时海洋环境因素的变化使滑翔机的运动参数存在一定的抖动，但二者的运动参数时间序列基本一致，说明了本文模型的正确性。

3.2 锯齿运动仿真

锯齿运动时设定控制参量{mb,rp} 的输入值为{0 .08 kg,0.02 m} ，初始运动参量{u0,v0,w0,p0,q0,r0,φ0,θ0,ψ0,x0,y0,z0} 均取值为0，每隔半个仿真周期将控制输入取反。令翼面形变量h在0～6 mm 范围内变化，仿真得到的不同翼面变形下滑翔机最终稳态速度、俯仰角、前进距离与最大下潜深度变化如图12所示，以0 mm和4 mm形变为例，两种形变下运动参量仿真结果随时间的变化如图13所示。

由图12 可知，在相同的控制输入下，锯齿运动模式时滑翔机下潜速度、俯仰角、下潜深度和水平位移总体呈下降趋势，说明翼形变对滑翔机锯齿运动的续航能力有所减弱。由图13 可知，滑翔机的速度、位移、姿态角等参数在控制量改变时存在短时间的波动，但滑翔机大部分时间处于稳定工作状态。4 mm 翼形变相对于无形变时滑翔机的最大前进距离减小了123.92 m（20.08%），最大的下潜深度减小了15.33 m（16.56%）。在一个仿真周期内，翼面的形变使滑翔机的续航能力和下潜深度降低，因此在实际航行过程中可主动调整滑翔机的控制参数去间接改善翼的形变，使翼面形变趋于更小。

3.3 螺旋运动仿真

螺旋运动时设定控制参量{γ,mb,rp} 的输入为{1 5°,0.05 kg,0.03 m} ，初始运动参量{u0,v0,w0,p0,q0,r0,φ0,θ0,ψ0,x0,y0,z0} 均取值为0，仿真时间设为1000 s，滑翔机达到稳定状态时的输出参量随形变的变化如图14 所示，以0 mm 和4 mm 翼面形变为例，两种形变下运动参量仿真结果随时间的变化如图15 所示。

由图14～15 的仿真结果可知，柔性可折叠翼水下滑翔机螺旋运动的回转半径、俯仰角、最大深度随翼面形变的增大而减小，说明翼面形变会降低滑翔机的续航能力，同时翼面形变使得滑翔机具有更小的螺旋回转半径，说明翼面形变使得滑翔机具有更好的机动性。在相同控制量输入下，4 mm 形变量相对于无形变量时滑翔机最大下潜深度减小了32.67 m（14.74%），回转半径减小了4.56 m（14.83%）。

4 结 论

本文对适用于载机空投和水下管式发射的柔性可折叠翼水下滑翔机水动力特性与滑翔运动性能进行了分析，通过比较柔性可折叠翼水下滑翔机不同翼面形变下的水动力特性、建立滑翔机的动力学模型并进行运动仿真，得出如下结论：

（1）翼面形变在0～6 mm 范围内，滑翔机俯仰力矩系数与形变量近似呈正弦变化关系，升力和阻力系数与形变量近似呈线性变化关系。

（2）动力学模型运动仿真结果与海试结果具有较好的一致性，验证了水动力数值计算和模型建立的准确性与可靠性。

（3）柔性可折叠翼水下滑翔机最大下潜深度和水平位移随翼面形变的增加而减小，翼面形变使螺旋运动模式下的回转半径减小。

上述结论可为柔性可折叠翼水下滑翔机控制决策提供依据。锯齿运动模式时，为保证滑翔机的续航能力和最大作业深度，可调整滑翔机的控制参数使翼面的形变趋于更小；螺旋运动模式时，可通过调整控制参数使翼面形变适当增加，以获得更小的螺旋回转半径进而提高滑翔机的机动性。