李振宇，陳質(zhì)二，俞建成，楊 慶，4

（1.中國科學(xué)院沈陽自動化研究所機(jī)器人學(xué)國家重點(diǎn)實(shí)驗(yàn)室，沈陽 110016；2.中國科學(xué)院機(jī)器人與智能制造創(chuàng)新研究院，沈陽 110169；3.中國科學(xué)院大學(xué)，北京 100049；4.東北大學(xué)機(jī)械工程與自動化學(xué)院，沈陽 110819）

0 引 言

水下滑翔機(jī)由于自身航速較低且續(xù)航時間有限[1]，其通常由母船送至目標(biāo)海域布放。在多數(shù)情況下觀測海區(qū)距離陸基較遠(yuǎn)，航渡時間較長（一般需要數(shù)天），尤其是在面對海上突發(fā)事件和海洋現(xiàn)象時滑翔機(jī)的響應(yīng)速度較慢，同時對于具有隱蔽布放需求的場合，如圖1所示，由于母船形體過大，布放活動易于暴露。因此通過載機(jī)空投和潛艇魚雷管發(fā)射來實(shí)現(xiàn)水下滑翔機(jī)的部署逐漸得到重視和研究。相比于傳統(tǒng)的剛性翼面水下滑翔機(jī)，柔性可折疊翼水下滑翔機(jī)的機(jī)翼在入水前折疊，入水后自動彈開，可以解決傳統(tǒng)滑翔機(jī)空投時剛性翼板存在的入水沖擊損壞[2]和潛艇水下發(fā)射管徑受限問題。

近年來，柔性可折疊翼較多地應(yīng)用在航空信息以及仿生領(lǐng)域，尤其是撲翼類微型飛行器方面。韓國建國大學(xué)的Muhammad 等[3]根據(jù)英國埃克塞特大學(xué)的Haas 等提出的四平面理論，利用形狀記憶合金絲作為驅(qū)動器，成功實(shí)現(xiàn)了可折疊翼的自動折疊展開動作；后來Truong 等[4]提出了基于四桿機(jī)構(gòu)的可折疊翼，該翼最大優(yōu)勢在于結(jié)構(gòu)簡單，缺點(diǎn)在于需要依靠外力推動連桿結(jié)構(gòu)才能實(shí)現(xiàn)折疊，翼在撲動中形態(tài)發(fā)生變化后無法及時修正；南京航空航天大學(xué)顧鑫[5]對柔性折疊翼飛行器飛行動力學(xué)進(jìn)行研究，分析了彈性變形和折疊角速度對縱向運(yùn)動參數(shù)的影響規(guī)律；斯坦福大學(xué)的Stowers等[6]仿照鳥類和蝙蝠的翅膀形態(tài)，設(shè)計(jì)了一種利用離心加速度使機(jī)翼從初始折疊狀態(tài)轉(zhuǎn)變?yōu)檎归_狀態(tài)的仿生可折疊翼；浙江大學(xué)李鐵風(fēng)等[7]研制出仿生軟體柔性機(jī)器魚，其由軟體人工肌肉驅(qū)動一對翼狀柔性胸鰭，通過節(jié)律性撲翅實(shí)現(xiàn)游動，翼翅選用耐高壓的軟體材料，可以承受萬米級別深海靜水壓力。

國內(nèi)外學(xué)者也對具有變翼功能的水下滑翔機(jī)開展了研究。大阪大學(xué)Arima 等[8-9]設(shè)計(jì)出可轉(zhuǎn)動翼板翻轉(zhuǎn)角度的“ALEX”滑翔機(jī)，其通過建立滑翔機(jī)的動力學(xué)模型并對不同翼板翻轉(zhuǎn)角度的滑翔機(jī)進(jìn)行運(yùn)動仿真，仿真結(jié)果表明，可轉(zhuǎn)動翼板角度的滑翔機(jī)具有更小的滑翔攻角和更好的運(yùn)動性能；天津大學(xué)楊志金等[10]提出一種具有變體功能的可變翼混合驅(qū)動水下滑翔機(jī)，通過改變機(jī)翼展弦比、后掠角和翻轉(zhuǎn)角可提高滑翔機(jī)在不同工況下的航行性能；西北工業(yè)大學(xué)田文龍等[11]對可變翼混合驅(qū)動水下滑翔機(jī)進(jìn)行研究，基于Newton-Euler方程建立了六自由度運(yùn)動數(shù)學(xué)模型，通過求解穩(wěn)態(tài)滑翔運(yùn)動平衡方程，定量分析了可控翼翻轉(zhuǎn)角對滑翔機(jī)滑翔性能的影響。上述多為機(jī)翼自由度數(shù)目增加的變翼研究，水下滑翔機(jī)通過載機(jī)空投和潛艇魚雷管發(fā)射時需要機(jī)翼能夠完全折疊。

本文以中國科學(xué)院沈陽自動化研究所自主研制的海翼1000mini水下滑翔機(jī)為柔性可折疊翼搭載基體，研究水下滑翔翼面展開后的水動力特性并建立其在流體中的動力學(xué)模型。根據(jù)水動力與運(yùn)動仿真實(shí)驗(yàn)的結(jié)果，分析翼形變與水下滑翔機(jī)水動力特性和運(yùn)動參數(shù)之間的關(guān)系，可為后續(xù)滑翔機(jī)控制算法設(shè)計(jì)和運(yùn)動優(yōu)化提供理論依據(jù)。

1 水動力參數(shù)辨識

1.1 翼面形位描述與形變量定義

柔性可折疊翼水下滑翔機(jī)機(jī)翼展開后柔性翼面受周圍流場的擾動會產(chǎn)生一定的形變，折疊翼的翼面支撐桿為不銹鋼材質(zhì)，相較于翼面可視其為剛性桿。根據(jù)剛?cè)狁詈隙囿w系統(tǒng)節(jié)點(diǎn)描述方法[12]，建立如圖2 所示的慣性坐標(biāo)系（I-i,j,k）和機(jī)翼坐標(biāo)系（S-x,y,z），翼面上某一點(diǎn)f0的形變?yōu)闄C(jī)翼坐標(biāo)系相對于慣性坐標(biāo)系的位置矢量和柔性翼面相對于機(jī)翼坐標(biāo)系變形矢量的疊加，對于翼面上某一點(diǎn)fi變形后對應(yīng)的Fi在慣性坐標(biāo)系下表示為

式中，m為Fi在慣性系中位置向量，K為機(jī)翼坐標(biāo)系S在慣性坐標(biāo)系I中對應(yīng)的位置向量，A為方向余弦矩陣，li為翼未發(fā)生形變時fi點(diǎn)在機(jī)翼坐標(biāo)系中的位置向量，hi為柔性翼的相對變形量。

在慣性系中的機(jī)翼坐標(biāo)系的位置向量M與方向余弦矩陣A均為定值，所以將形變描述簡化至機(jī)翼坐標(biāo)系中，下文中均選定翼面外邊緣弧頂至翼面坐標(biāo)系y軸的標(biāo)量距離h=‖h‖作為翼面的形變量。

1.2 水動力修正計(jì)算公式引入

水下滑翔機(jī)航行時所受到的水動力包括滑翔機(jī)的升力L、阻力D與俯仰力矩M，當(dāng)滑翔機(jī)以小攻角航行，不考慮翼形變時滑翔機(jī)的升力、阻力和力矩滿足如下關(guān)系式[13]：

式中，ρ為海水密度，A為橫截面積，V為航行速度，α為航行攻角，CD為阻力系數(shù)，CL為升力系數(shù)，CM為俯仰力矩系數(shù)。KD0、KL0、KM0為與滑翔機(jī)直航特性有關(guān)的水動力系數(shù)，KD、KL、KM為與滑翔機(jī)攻角特性有關(guān)的水動力系數(shù)。

翼形變h會對滑翔機(jī)升力、阻力和俯仰力矩產(chǎn)生影響，可將式（2）中水動力系數(shù)寫為與翼形變有關(guān)的變量，因此將柔性可折疊翼水下滑翔機(jī)水動力計(jì)算公式修正如下：

式 中，KWD0(h)、KWD(h)、KWL0(h)、KWL(h)、KWM0(h)、KWM(h)為修正后的與滑翔機(jī)翼面形變有關(guān)的水動力系數(shù)。

1.3 柔性可折疊翼結(jié)構(gòu)設(shè)計(jì)與工作原理

如圖3 所示，柔性可折疊水下滑翔機(jī)在其艉部密封鋼圈和機(jī)身回轉(zhuǎn)體上安裝一對可折疊機(jī)翼。其中，翼面支撐桿底端轉(zhuǎn)動塊和基體鉸接在一起，鉸接軸外圍并排安裝三個扭簧，扭簧驅(qū)動轉(zhuǎn)動塊繞鉸接軸旋轉(zhuǎn)。支撐桿底端可轉(zhuǎn)動體和基體之間設(shè)計(jì)一對棘輪自鎖機(jī)構(gòu)，使開合后的翼處于目標(biāo)位置且不會受外界影響產(chǎn)生復(fù)位和抖動?；铏C(jī)入水前用水溶性膠帶捆綁使其處于收合狀態(tài)的機(jī)翼，滑翔機(jī)入水后水溶性膠帶遇水溶解斷裂，折疊翼迅速自動展開，水溶膠帶斷裂瞬間翼面展開過程如圖4所示。

1.4 水動力仿真實(shí)驗(yàn)設(shè)計(jì)

本文通過STAR-CCM+水動力仿真軟件構(gòu)建仿真模型，計(jì)算流域是長為10l，寬和高均為3.5l的長方體（l為滑翔機(jī)的長度），如圖5所示?；铏C(jī)主軸線與流域的長軸重合，流域速度進(jìn)口端為長方體的前端面，與滑翔機(jī)艏部距離為4l，流域壓力出口端為長方體的后端面，與滑翔機(jī)艉部距離為5l。計(jì)算水域采用多面體網(wǎng)格進(jìn)行網(wǎng)格劃分，其CFD 分析模型的網(wǎng)格離散方式如圖6 所示，為提高計(jì)算精度，采用棱柱層網(wǎng)格對水下滑翔機(jī)近壁面湍流進(jìn)行捕捉，其中棱柱層網(wǎng)格近壁厚度設(shè)置為1 mm，層數(shù)設(shè)為8 層，網(wǎng)格增長率為1.3。仿真模型中加入面控制對滑翔機(jī)的機(jī)體和翼的表面進(jìn)行了加密。建立四組不同網(wǎng)格尺寸的計(jì)算模型開展網(wǎng)格無關(guān)性驗(yàn)證，四組網(wǎng)格尺寸下網(wǎng)格數(shù)量分別為80萬、120萬、150萬、200萬，仿真計(jì)算均在1.0 m/s航速和4°攻角下進(jìn)行，不同網(wǎng)格數(shù)量下阻力的仿真計(jì)算結(jié)果如表1所示。通過對比分析可知，當(dāng)網(wǎng)格數(shù)達(dá)到150萬時滑翔機(jī)的阻力基本上不再發(fā)生變化，綜合考慮計(jì)算精度和效率，網(wǎng)格數(shù)量選取為150萬。

表1 不同網(wǎng)格數(shù)量阻力計(jì)算結(jié)果Tab.1 Computational results of lift-drag ratio of different mesh numbers

已知流體運(yùn)動粘性系數(shù)υ= 1.003 74×10-6/(m2·s-1)，設(shè)定滑翔機(jī)最大航速為1 m/s，機(jī)身總長度為1.08 m，求得對應(yīng)的雷諾系數(shù)為

所以計(jì)算過程應(yīng)選擇湍流模型，由于Re在106數(shù)量級，本文在Star CCM+仿真軟件中選擇低雷諾數(shù)k-ε湍流模型進(jìn)行流體仿真計(jì)算。

1.5 水動力計(jì)算公式的確定

為確定不同翼形變量時柔性可折疊翼水下滑翔機(jī)水動力特性變化規(guī)律，將不同形變后的翼面簡化為剛性翼處理，對0～6 mm 形變時滑翔機(jī)的水動力特性開展仿真實(shí)驗(yàn)，以0.5 m/s 航速下水動力實(shí)驗(yàn)結(jié)果為例，水動力數(shù)值計(jì)算結(jié)果如圖7所示。

由圖可知，翼面存在形變時滑翔機(jī)升力、阻力、俯仰力矩與攻角分別滿足線性、二次型、線性關(guān)系，符合海翼系列水下滑翔機(jī)不同攻角下的水動力變化規(guī)律[14]。其中，升力和阻力同時隨翼面形變的增加而減小，俯仰力矩隨翼形變量增加先變大后變小。翼面形變的增加使0°～8°攻角下對應(yīng)的升阻比逐漸減小，當(dāng)攻角大于10°時，不同翼面形變下的升阻比基本保持不變。根據(jù)圖7 中得到水動力實(shí)驗(yàn)仿真結(jié)果分析式（3）中升力、阻力、俯仰力矩對應(yīng)的水動力系數(shù)KWL0、KWL、KWD0、KWD、KWM0、KWM與形變量之間的關(guān)系，各水動力系數(shù)擬合計(jì)算結(jié)果如圖8 所示，可以看出KWM0、KWM與形變量近似呈正弦關(guān)系，KWD0、KWD、KWL、KWL0與形變量近似呈線性關(guān)系，則柔性可折疊翼水下滑翔機(jī)水動力計(jì)算滿足下式：

用( · )表示D、L、M，則( ·)a0-d0、( ·)a-d分別表示水動力計(jì)算修正公式中直航特性與攻角特性水動力系數(shù)，將擬合后的結(jié)果代入上式，得到修正后的柔性可折疊翼滑翔機(jī)水動力計(jì)算公式為

2 動力學(xué)建模

2.1 坐標(biāo)系的建立與坐標(biāo)變換

如圖9所示，建立慣性坐標(biāo)系E0-ijk，向量b=[x,y,z]T和向量S=[φ,θ,ψ]T表示滑翔機(jī)的浮心在慣性坐標(biāo)系中的位置和姿態(tài)。再以滑翔機(jī)浮心為坐標(biāo)原點(diǎn)建立機(jī)體坐標(biāo)系e0-e1e2e3，分別用V=[u,v,w]T和Ω=[p,q,r]T表示滑翔機(jī)在機(jī)體坐標(biāo)系下的速度與角速度。

速度坐標(biāo)系π0-π1π2π3原點(diǎn)與機(jī)體坐標(biāo)系原點(diǎn)重合，令機(jī)體坐標(biāo)系先繞e2軸轉(zhuǎn)動α，使得e3與π3重合，然后再將旋轉(zhuǎn)之后的坐標(biāo)系繞π3轉(zhuǎn)動β，此時e1與e2分別與速度坐標(biāo)系坐標(biāo)軸π1和π2重合。所以從速度坐標(biāo)系到機(jī)體坐標(biāo)系之間的轉(zhuǎn)換矩陣可表示為

根據(jù)慣性坐標(biāo)系和機(jī)體坐標(biāo)系之間的轉(zhuǎn)換關(guān)系[15]可知

2.2 動力學(xué)模型

水下滑翔機(jī)是由殼體、可移動質(zhì)量塊、姿態(tài)調(diào)節(jié)系統(tǒng)等組成的多剛體系統(tǒng)[14]。對于圖10所示的柔性可折疊翼滑翔機(jī)，將其整體分解為三個質(zhì)點(diǎn)：殼體質(zhì)點(diǎn)（質(zhì)量ms，質(zhì)心rs）、可移動電池塊質(zhì)點(diǎn)（質(zhì)量mp，質(zhì)心rp）、浮力調(diào)節(jié)質(zhì)點(diǎn)（質(zhì)量mb，質(zhì)心rb）。由上述劃分可知系統(tǒng)的總質(zhì)量為

已知電池塊質(zhì)心偏心距為Rp，當(dāng)其轉(zhuǎn)動角度為γ，沿e1軸移動距離為rp時，其質(zhì)心rp可表示為

令p和P分別表示滑翔機(jī)在慣性坐標(biāo)系下和機(jī)體坐標(biāo)系下的動量，π和Π表示慣性坐標(biāo)系和機(jī)體坐標(biāo)系下的角動量，md表示排水量。參照牛頓-歐拉方程，滑翔機(jī)系統(tǒng)的動量和受力滿足如下關(guān)系[13]：

式中，fext和τext分別為外界環(huán)境對滑翔機(jī)產(chǎn)生的力和力矩，lsE、lpE、lbE分別為圖10 中慣性坐標(biāo)系下殼體質(zhì)點(diǎn)、電池塊質(zhì)點(diǎn)、浮力質(zhì)點(diǎn)對應(yīng)的位置矢量，g表示重力加速度。

定義運(yùn)算符·為向量“·”的叉乘，根據(jù)慣性坐標(biāo)系和體坐標(biāo)之間坐標(biāo)變換關(guān)系進(jìn)一步推導(dǎo)可得柔性可折疊翼滑翔機(jī)廣義力和力矩在體系中的表達(dá)式為

其中，F(xiàn)W=Rfext，TW=Rτext，二者為柔性可折疊翼滑翔機(jī)在流體中受到的粘性水動力項(xiàng)，根據(jù)式（5）可知，在體坐標(biāo)系表達(dá)式為

式中，V=‖V‖，表示機(jī)體坐標(biāo)系中滑翔機(jī)相對于海流的速度，KMR、Kp、Kβ、Kq、KMY、Kr為滑翔機(jī)自身的水動力系數(shù)。

定義M為系統(tǒng)的廣義質(zhì)量矩陣，其由附加質(zhì)量矩陣Mf=diag([Mf1Mf2Mf3])、附加轉(zhuǎn)動慣量矩陣IA=diag([If1If2If3])和耦合項(xiàng)組成。將柔性可折疊翼滑翔機(jī)機(jī)體坐標(biāo)系下的運(yùn)動特征用矢量u=[V,Ω]表述，動量特征用η=[P Π]表示，則η和u之間的關(guān)系可表示為

上式對時間進(jìn)行求導(dǎo)，可得

根據(jù)上式等可推出u?的表達(dá)式為

再將上式進(jìn)一步展開得到柔性可折疊翼水下滑翔機(jī)運(yùn)動方程為

3 仿真與實(shí)驗(yàn)分析

根據(jù)滑翔機(jī)動力學(xué)模型編寫仿真程序并使用四階龍格庫塔法對式（20）微分方程組進(jìn)行求解[17-18]，得到滑翔機(jī)滑翔運(yùn)動時各運(yùn)動參量隨時間的變化。仿真分析時滑翔機(jī)的水動力系數(shù)和設(shè)計(jì)參數(shù)如表2～3所示。

表2 柔性可折疊水下滑翔機(jī)水動力系數(shù)Tab.2 Hydrodynamic coefficients of the flexible foldable-wings glider

表3 柔性可折疊翼水下滑翔機(jī)設(shè)計(jì)參數(shù)Tab.3 Design parameters of the flexible foldable-wings glider

3.1 仿真與海試比較分析

為驗(yàn)證本文計(jì)算和仿真的正確性，在輸入相同控制量條件下，將柔性可折疊翼水下滑翔機(jī)0形變仿真結(jié)果與海翼1000mini水下滑翔機(jī)2020年12月南海海試實(shí)驗(yàn)數(shù)據(jù)進(jìn)行比對分析，滑翔機(jī)的深度和俯仰角的變化如圖11 所示，海試時海洋環(huán)境因素的變化使滑翔機(jī)的運(yùn)動參數(shù)存在一定的抖動，但二者的運(yùn)動參數(shù)時間序列基本一致，說明了本文模型的正確性。

3.2 鋸齒運(yùn)動仿真

鋸齒運(yùn)動時設(shè)定控制參量{mb,rp} 的輸入值為{0 .08 kg,0.02 m} ，初始運(yùn)動參量{u0,v0,w0,p0,q0,r0,φ0,θ0,ψ0,x0,y0,z0} 均取值為0，每隔半個仿真周期將控制輸入取反。令翼面形變量h在0～6 mm 范圍內(nèi)變化，仿真得到的不同翼面變形下滑翔機(jī)最終穩(wěn)態(tài)速度、俯仰角、前進(jìn)距離與最大下潛深度變化如圖12所示，以0 mm和4 mm形變?yōu)槔瑑煞N形變下運(yùn)動參量仿真結(jié)果隨時間的變化如圖13所示。

由圖12 可知，在相同的控制輸入下，鋸齒運(yùn)動模式時滑翔機(jī)下潛速度、俯仰角、下潛深度和水平位移總體呈下降趨勢，說明翼形變對滑翔機(jī)鋸齒運(yùn)動的續(xù)航能力有所減弱。由圖13 可知，滑翔機(jī)的速度、位移、姿態(tài)角等參數(shù)在控制量改變時存在短時間的波動，但滑翔機(jī)大部分時間處于穩(wěn)定工作狀態(tài)。4 mm 翼形變相對于無形變時滑翔機(jī)的最大前進(jìn)距離減小了123.92 m（20.08%），最大的下潛深度減小了15.33 m（16.56%）。在一個仿真周期內(nèi)，翼面的形變使滑翔機(jī)的續(xù)航能力和下潛深度降低，因此在實(shí)際航行過程中可主動調(diào)整滑翔機(jī)的控制參數(shù)去間接改善翼的形變，使翼面形變趨于更小。

3.3 螺旋運(yùn)動仿真

螺旋運(yùn)動時設(shè)定控制參量{γ,mb,rp} 的輸入為{1 5°,0.05 kg,0.03 m} ，初始運(yùn)動參量{u0,v0,w0,p0,q0,r0,φ0,θ0,ψ0,x0,y0,z0} 均取值為0，仿真時間設(shè)為1000 s，滑翔機(jī)達(dá)到穩(wěn)定狀態(tài)時的輸出參量隨形變的變化如圖14 所示，以0 mm 和4 mm 翼面形變?yōu)槔?，兩種形變下運(yùn)動參量仿真結(jié)果隨時間的變化如圖15 所示。

由圖14～15 的仿真結(jié)果可知，柔性可折疊翼水下滑翔機(jī)螺旋運(yùn)動的回轉(zhuǎn)半徑、俯仰角、最大深度隨翼面形變的增大而減小，說明翼面形變會降低滑翔機(jī)的續(xù)航能力，同時翼面形變使得滑翔機(jī)具有更小的螺旋回轉(zhuǎn)半徑，說明翼面形變使得滑翔機(jī)具有更好的機(jī)動性。在相同控制量輸入下，4 mm 形變量相對于無形變量時滑翔機(jī)最大下潛深度減小了32.67 m（14.74%），回轉(zhuǎn)半徑減小了4.56 m（14.83%）。

4 結(jié) 論

本文對適用于載機(jī)空投和水下管式發(fā)射的柔性可折疊翼水下滑翔機(jī)水動力特性與滑翔運(yùn)動性能進(jìn)行了分析，通過比較柔性可折疊翼水下滑翔機(jī)不同翼面形變下的水動力特性、建立滑翔機(jī)的動力學(xué)模型并進(jìn)行運(yùn)動仿真，得出如下結(jié)論：

（1）翼面形變在0～6 mm 范圍內(nèi)，滑翔機(jī)俯仰力矩系數(shù)與形變量近似呈正弦變化關(guān)系，升力和阻力系數(shù)與形變量近似呈線性變化關(guān)系。

（2）動力學(xué)模型運(yùn)動仿真結(jié)果與海試結(jié)果具有較好的一致性，驗(yàn)證了水動力數(shù)值計(jì)算和模型建立的準(zhǔn)確性與可靠性。

（3）柔性可折疊翼水下滑翔機(jī)最大下潛深度和水平位移隨翼面形變的增加而減小，翼面形變使螺旋運(yùn)動模式下的回轉(zhuǎn)半徑減小。

上述結(jié)論可為柔性可折疊翼水下滑翔機(jī)控制決策提供依據(jù)。鋸齒運(yùn)動模式時，為保證滑翔機(jī)的續(xù)航能力和最大作業(yè)深度，可調(diào)整滑翔機(jī)的控制參數(shù)使翼面的形變趨于更小；螺旋運(yùn)動模式時，可通過調(diào)整控制參數(shù)使翼面形變適當(dāng)增加，以獲得更小的螺旋回轉(zhuǎn)半徑進(jìn)而提高滑翔機(jī)的機(jī)動性。