飛機(jī)圓周盤旋時(shí)拖纜的動(dòng)力學(xué)建模與分析

朱延娟， 羅 丹

(同濟(jì)大學(xué) 航空航天與力學(xué)學(xué)院，上海 200092)

航空拖曳系統(tǒng)主要是由飛機(jī)、拖纜和拖體三者組成的多體約束系統(tǒng)，其中拖纜是飛機(jī)與拖體之間的唯一連接介質(zhì)，主要起能量傳遞和電信號(hào)傳輸?shù)淖饔肹1].航空拖曳系統(tǒng)有著許多實(shí)際的工程應(yīng)用，如防空武器系統(tǒng)鑒定、軍事作戰(zhàn)目標(biāo)模擬和通信系統(tǒng)等領(lǐng)域.飛機(jī)通信應(yīng)用是目前各國(guó)重點(diǎn)研究的關(guān)鍵技術(shù)問(wèn)題，拖纜的穩(wěn)態(tài)構(gòu)型和垂直度等是檢驗(yàn)系統(tǒng)是否正常有效工作的重要性能指標(biāo)，是設(shè)計(jì)中必須重點(diǎn)研究的主要問(wèn)題.垂直度是拖纜在鉛垂線上的投影長(zhǎng)度與總長(zhǎng)度之比，為確保系統(tǒng)能夠正常有效工作，垂直度需不小于70%.若要保證垂直度在要求范圍內(nèi)，飛機(jī)必須在特定狀態(tài)下進(jìn)行穩(wěn)定盤旋飛行[2].

1958年，Caughey[3]分析了拖索系統(tǒng)在真空中的運(yùn)動(dòng)，提出了動(dòng)態(tài)穩(wěn)定性的近似求解方法，但忽略了空氣阻尼對(duì)拖索系統(tǒng)的影響.Lemon和Fraser[4]考慮空氣阻力的情況，假設(shè)拖索是不可伸長(zhǎng)的，發(fā)現(xiàn)在低阻尼和高旋轉(zhuǎn)速度下拖索系統(tǒng)會(huì)存在不穩(wěn)定性.Skop和Choo[5]把末端連接物體假設(shè)為球體，研究柔性不可伸長(zhǎng)拖索的平衡條件與構(gòu)型，該研究?jī)H考慮了法向空氣阻力.Huang[6]在圓柱坐標(biāo)系下，根據(jù)牛頓第二定律推導(dǎo)拖索的動(dòng)力學(xué)方程，確定了平衡狀態(tài)的拖索構(gòu)型，此方法對(duì)初始條件有著較高的要求.Zhu和Rahn[7]建立了拖索系統(tǒng)的動(dòng)力學(xué)模型，分析系統(tǒng)的動(dòng)態(tài)響應(yīng)與穩(wěn)定性，得出在小阻尼、大質(zhì)量末端連接物和高旋轉(zhuǎn)速度下存在穩(wěn)態(tài)多解.鄭小洪等[8]應(yīng)用牛頓定律建立了機(jī)載甚低頻拖曳天線的動(dòng)力學(xué)模型，運(yùn)用Newton-Raphson迭代法求解天線穩(wěn)定垂直構(gòu)型，分析了空氣系數(shù)和末端質(zhì)量對(duì)天線穩(wěn)定構(gòu)型的影響.此外，Choo和Casarella[9]對(duì)拖索動(dòng)態(tài)仿真方法進(jìn)行比較分析，認(rèn)為有限元是較為通用的動(dòng)態(tài)仿真法.Russell和Anderson[10]運(yùn)用有限元建立了拖索系統(tǒng)模型，把末端物體簡(jiǎn)化為質(zhì)點(diǎn)，定性研究了拖索系統(tǒng)的穩(wěn)態(tài)平衡.James和Louis[11]在航空拖索上建立局部坐標(biāo)系，應(yīng)用有限元將拖纜離散成數(shù)分段，針對(duì)任一拖纜段建立平衡方程，運(yùn)用中心差分法對(duì)其進(jìn)行求解與分析.Paul和Pavel[12]應(yīng)用有限元對(duì)拖索進(jìn)行了建模研究，把剛性段作為無(wú)質(zhì)量的二力桿，段與段之間用球鉸連接，每段拖索的質(zhì)量和受力都集中在節(jié)點(diǎn)上.拖索離散模型的研究能夠有效地處理強(qiáng)非線性問(wèn)題，并在一定程度上反映拖索的實(shí)際工作情況.目前研究拖索的相關(guān)文獻(xiàn)較多局限于分析拖索的受力情況與方程求解.

在上述研究基礎(chǔ)上，針對(duì)飛機(jī)的圓周盤旋運(yùn)動(dòng)，采用多剛體法建立了拖纜離散模型.通過(guò)對(duì)拖纜模型的數(shù)學(xué)分析，構(gòu)建了一種便于仿真計(jì)算的拖纜多剛體動(dòng)力學(xué)模型.并對(duì)拖纜連續(xù)模型應(yīng)用微元法進(jìn)行理論驗(yàn)證分析，結(jié)果表明所建的多剛體動(dòng)力學(xué)模型是正確可信的.

1 拖纜的動(dòng)力學(xué)建模

1.1 相關(guān)坐標(biāo)系

對(duì)拖纜模型進(jìn)行力學(xué)分析，需構(gòu)建空間坐標(biāo)系.空間坐標(biāo)系由慣性坐標(biāo)系OXYZ和旋轉(zhuǎn)坐標(biāo)系Oxyz構(gòu)成，兩者具有相同原點(diǎn)，旋轉(zhuǎn)坐標(biāo)系繞縱軸以飛機(jī)旋轉(zhuǎn)角速度旋轉(zhuǎn).將拖纜離散為任意n段，運(yùn)用有限元中的剛體法建立拖纜的多剛體模型.任意段拖纜有著獨(dú)立的旋轉(zhuǎn)坐標(biāo)系R(x,y,z)和描述旋轉(zhuǎn)位置方向的歐拉角γ(ψ,θ,φ).設(shè)旋轉(zhuǎn)坐標(biāo)系和歐拉角作為廣義坐標(biāo)，即：

qi=[xi,yi,zi,ψi,θi,φi]T,i=1,2,…,n

1.2 受力分析

式中:E為楊氏模量；A0為無(wú)應(yīng)力作用下拖纜的橫截面積；Lsj、lj分別為第j段拖纜的無(wú)應(yīng)力長(zhǎng)度和應(yīng)力長(zhǎng)度.

重力方程：

Gj=-mjgK

式中:mj為第j段拖纜質(zhì)量；g為重力加速度；K為沿OZ軸方向的單位矢量.

圖1 拖纜段的受力示意圖

式中:CDj、CLj分別為第j段拖纜的空氣阻力系數(shù)、升力系數(shù);ρj為空氣密度;dj為拖纜段直徑;vj為拖纜與空氣的相對(duì)運(yùn)動(dòng)速度;eDj、eLj分別為空氣阻力和升力方向的單位矢量.

由于拖纜做三維旋轉(zhuǎn)運(yùn)動(dòng)時(shí)阻力方向不斷改變，在動(dòng)力學(xué)仿真計(jì)算中較難對(duì)其進(jìn)行準(zhǔn)確的設(shè)置與計(jì)算.因此，在拖纜段旋轉(zhuǎn)坐標(biāo)系上將總空氣阻力分解為沿慣性坐標(biāo)系X、Y、Z軸方向的阻力分量，如圖2所示.

拖纜力學(xué)模型的合理構(gòu)建，能快速、準(zhǔn)確地求解出運(yùn)動(dòng)狀態(tài).在三維流場(chǎng)中，第j段拖纜的空氣阻力分量Fxj、Fyj、Fzj的方程式[14]如下：

式中：Vxj、Vyj、Vzj分別為第j段拖纜沿慣性坐標(biāo)系x、y、z方向的流場(chǎng)速度；π為圓周率；cf、cp分別為空氣摩擦系數(shù)和壓力系數(shù).

圖2 拖纜段的力學(xué)模型

1.3 運(yùn)動(dòng)方程

拖纜的運(yùn)動(dòng)規(guī)律由運(yùn)動(dòng)方程求解得出.拖纜系統(tǒng)模型的自由度為3n(n為拖纜的段數(shù))，變量為廣義坐標(biāo)，則拖纜的拉格朗日運(yùn)動(dòng)方程為

第j段拖纜的慣性位置：

式中：C(ψj,θj,φj)為坐標(biāo)轉(zhuǎn)換矩陣；Xj、Yj、Zj分別為第j段拖纜的慣性坐標(biāo)位移；Lsj為第i段拖纜的無(wú)應(yīng)力長(zhǎng)度.

總動(dòng)能T是所有拖纜段動(dòng)能的總和，則：

式中：mj，Vj分別為第j段拖纜的質(zhì)量和速度.

i=1,2,…,3n

式中：FXj、FYj、FZj是第j段拖纜的慣性作用力.

2 拖纜運(yùn)動(dòng)的理論驗(yàn)證分析

多剛體法建立拖纜離散模型，可以省去拖纜段間約束力的繁瑣計(jì)算，動(dòng)力學(xué)模型的合理構(gòu)建，能夠快速準(zhǔn)確獲得拖纜的運(yùn)動(dòng)狀態(tài).由于拖纜連續(xù)模型更接近于實(shí)際形狀，因此對(duì)拖纜連續(xù)模型應(yīng)用微元法進(jìn)行理論分析，并與拖纜離散系統(tǒng)的穩(wěn)態(tài)構(gòu)型結(jié)果進(jìn)行對(duì)比.

2.1 連續(xù)模型的構(gòu)建

采用有限元法是把拖纜離散化成一系列的拖纜段模型，而微元法是在拖纜任意一點(diǎn)上取細(xì)小微元進(jìn)行分析.設(shè)拖纜自然無(wú)伸展?fàn)顟B(tài)、穩(wěn)定狀態(tài)和動(dòng)態(tài)平衡狀態(tài)下任意位置E的弧坐標(biāo)分別為S0、S1、S2.任意E點(diǎn)的位移表示為

U(S0，T)=R2(S0，T)-R1(S0)

(1)

式中:U(S0,T)為E點(diǎn)從穩(wěn)定狀態(tài)到動(dòng)態(tài)平衡T時(shí)刻的位移；R2(S0,T)為動(dòng)態(tài)平衡中O點(diǎn)到E點(diǎn)在T時(shí)刻的距離；R1(S0)為穩(wěn)定狀態(tài)中O點(diǎn)到E點(diǎn)的距離.

在拖纜任意E點(diǎn)取微元，對(duì)其進(jìn)行受力分析如圖3所示，根據(jù)牛頓定律可得拖纜的運(yùn)動(dòng)方程：

(2)

式中:P為拖纜所受的張力;Fw為拖纜重力;FI為拖纜的慣性力;FD為法向空氣阻力.

我們?yōu)槭裁磿?huì)主張社會(huì)主義市場(chǎng)經(jīng)濟(jì)?這不僅是我們觀察世界的結(jié)果，也是回歸本民族歷史的思考。經(jīng)濟(jì)學(xué)理論告訴我們，好的制度應(yīng)該是“內(nèi)生”的，中國(guó)在社會(huì)主義市場(chǎng)經(jīng)濟(jì)體制下表現(xiàn)良好，是因?yàn)檫@種制度選擇符合我們自己的歷史。按照市場(chǎng)經(jīng)濟(jì)等于資本主義的標(biāo)準(zhǔn)，很多人不認(rèn)為中國(guó)古代是市場(chǎng)經(jīng)濟(jì)國(guó)家，我們?cè)诖饲也徽撨@種看法的偏頗性，我們僅僅強(qiáng)調(diào)，市場(chǎng)經(jīng)濟(jì)存在于古代中國(guó)，且為古代中國(guó)的發(fā)展做出了巨大貢獻(xiàn)。我們認(rèn)為，市場(chǎng)經(jīng)濟(jì)作為一種自發(fā)秩序，在中國(guó)和歐洲都有遙遠(yuǎn)的歷史。只要人們生產(chǎn)有了剩余需要交換，市場(chǎng)和貿(mào)易就會(huì)自然而然產(chǎn)生。當(dāng)然古代的市場(chǎng)規(guī)模與所包含的內(nèi)容與現(xiàn)代社會(huì)是不能比擬的，但從促進(jìn)資源配置的功能看卻是一致的。

圖3 拖纜微元受力示意圖

由于拖纜的切向阻力部分很小，可以忽略不計(jì)，因此FD只考慮法向空氣阻力[8].則：

(3)

式中：P為拖纜張力的數(shù)值大??；ρ1為拖纜單位長(zhǎng)度密度；ρ2為空氣密度；aE為任意E點(diǎn)的加速度；Dn為法向空氣阻力系數(shù)；Vn為拖纜相對(duì)于空氣的法向速度;A0為無(wú)應(yīng)力作用下拖纜的橫截面積;k為Oz軸的單位矢量.

拖纜末端連接著拖體，由于拖體的作用，因此存在固有邊界條件：

P+FDw+FDI+FDH=0

(4)

其中，拖體的重力FDw、慣性力FDI和空氣阻力FDH計(jì)算式為

(5)

式中:W=Mg、M、aD、D、VD分別是拖體的重力、質(zhì)量、加速度、空氣阻力系數(shù)和速度.

2.2 穩(wěn)態(tài)動(dòng)力學(xué)模型

在穩(wěn)態(tài)時(shí)，拖纜位移U=0，u=U/a,其中,a為飛機(jī)盤旋半徑.即u=0，將U=0、式(1)和式(3)代入方程式(2)，通過(guò)量綱一化得：

ρk+ω*2[k×(k×r)+dn|vn|vn]

(6)

其中，s=S0/a，p=P/ρ2A0ag，r=xi+yi+zk=R1/a，ρ=ρ1/ρ2，ω*2=aω2/g，t=ωT，dn=Dn/ρ2A0，vn=dr/ds×((k×r)×dr/ds).

在穩(wěn)態(tài)時(shí)，拖體位移Ud=U(0,T)=0，ud=Ud/a,即ud=0，將Ud=0和式(5)代入方程式(4)，量綱一化得：

wk+mω*2k×(k×r)+dω*2|k×r|(k×r)

(7)

其中，pd=p(0),即pd為s=0時(shí)拖纜的張力，w=W/ρ2aA0g，m=M/ρ2aA0，d=D/ρ2A0.

由方程式(6)、(7)得到：

上述量綱一模型的邊界條件如下：

在s=0處，

x=Xd/a=xd，y=z=0

式中:Xd為假定拖體在旋轉(zhuǎn)坐標(biāo)系下的初始Ox軸坐標(biāo);xd為計(jì)算的假定初值.

在s=L/a=l處，

x2+y2=1，u=U/a

3 拖纜的應(yīng)用分析

實(shí)際應(yīng)用中，飛機(jī)圓周盤旋處于甚低頻狀態(tài)，拖纜的活動(dòng)區(qū)域?qū)儆趯?duì)流層.考慮研究需求和結(jié)合實(shí)際情況，假設(shè)飛機(jī)圓周速度v=80 m·s-1；拖纜長(zhǎng)度l=5 000 m；拖纜密度ρ1=7 850 kg·m-3；拖纜直徑d=2 mm；拖體質(zhì)量M=50 kg；空氣密度ρ2=1.225 kg·m-3；空氣摩擦系數(shù)cf=0.02；空氣壓力系數(shù)cp=1.2；拖纜離散化數(shù)n=20段.

在拖纜系統(tǒng)工作時(shí)，飛機(jī)以特定的狀態(tài)進(jìn)行穩(wěn)定盤旋飛行[15]，飛機(jī)速度保持不變，拖纜在一定范圍內(nèi)擺動(dòng)，拖體保持在中心位置，即悠悠效應(yīng).

當(dāng)系統(tǒng)穩(wěn)定有效工作時(shí)，拖纜進(jìn)行動(dòng)態(tài)穩(wěn)定運(yùn)動(dòng)，穩(wěn)定構(gòu)型的仿真計(jì)算結(jié)果是呈螺旋形的垂直線，其俯視圖如圖4所示.從圖中可以看出隨著飛機(jī)的運(yùn)動(dòng)，拖纜在水平面內(nèi)進(jìn)行穩(wěn)定的盤旋運(yùn)動(dòng).基于打靶法應(yīng)用四階龍格庫(kù)塔迭代法，對(duì)拖纜連續(xù)模型的穩(wěn)態(tài)動(dòng)力學(xué)模型進(jìn)行數(shù)值求解，在相同參數(shù)下得到量綱一化構(gòu)型的俯視圖如圖5所示.可以發(fā)現(xiàn)兩種方法的構(gòu)型結(jié)果一致，說(shuō)明所構(gòu)建的多剛體動(dòng)力學(xué)模型與仿真計(jì)算是正確可信的.

圖4 拖纜系統(tǒng)穩(wěn)定構(gòu)型俯視圖

基于飛機(jī)的實(shí)際飛行情況，對(duì)飛機(jī)先直線飛行再轉(zhuǎn)入盤旋軌道運(yùn)動(dòng)進(jìn)行仿真計(jì)算，在此過(guò)程中不同參數(shù)下拖纜垂直長(zhǎng)度的變化情況如圖6～10所示.從圖中可知，飛機(jī)狀態(tài)轉(zhuǎn)換時(shí)，拖纜經(jīng)歷了上下大幅擺動(dòng)的過(guò)程，垂直長(zhǎng)度的變化過(guò)程類似于正余弦波動(dòng)圖，而后逐漸進(jìn)入工作的運(yùn)動(dòng)狀態(tài).結(jié)果顯示仿真600 s過(guò)后，拖纜進(jìn)行動(dòng)態(tài)穩(wěn)定運(yùn)動(dòng)，其垂直長(zhǎng)度將趨于保持不變.

不同飛機(jī)速度下拖纜垂直長(zhǎng)度的變化情況如圖6所示，當(dāng)60 m·s-190 m·s-1時(shí)，系統(tǒng)難以正常工作.不同拖纜直徑下拖纜垂直長(zhǎng)度的變化情況如圖7所示，當(dāng)2 mm

圖5 量綱一化拖纜構(gòu)型俯視圖

Fig.5Planeviewoftowedcableconfigurationindimensionlessquantitymethod

圖6 不同飛機(jī)速度下拖纜垂直長(zhǎng)度的變化情況

Fig.6Variationofverticallengthoftowedcableatdifferentspeeds

圖7 不同拖纜直徑下拖纜垂直長(zhǎng)度的變化情況

Fig.7Variationofverticallengthoftowedcableindifferenttowedcablediameter

圖8 不同拖體質(zhì)量下拖纜垂直長(zhǎng)度的變化情況

Fig.8Variationofverticallengthoftowedcableatdifferentmassoftowedbody

圖9 不同空氣密度下拖纜垂直長(zhǎng)度的變化情況

Fig.9Variationofverticallengthoftowedcableatdifferentairdensity

不同參數(shù)下拖纜垂直度大小如表1所示.可知增加末端拖體質(zhì)量類似于減少空氣阻尼的效果.除飛機(jī)速度外，其他參數(shù)在一定范圍內(nèi)變化，雖垂直度不同，但都處于系統(tǒng)正常工作范圍內(nèi).飛機(jī)速度臨界值在80～90 m·s-1范圍內(nèi)，則應(yīng)用二分法得出v≤83 m·s-1時(shí)滿足通信要求.因此，系統(tǒng)要能正常有效通信工作，必須選擇合適的飛機(jī)速度.

圖10 不同空氣系數(shù)下拖纜垂直長(zhǎng)度的變化情況

Fig.10Variationofverticallengthoftowedcableatdifferentaircoefficient

表1 不同參數(shù)下拖纜垂直度

4 結(jié)論

(1)針對(duì)飛機(jī)的圓周盤旋運(yùn)動(dòng)，構(gòu)建了拖纜的多剛體動(dòng)力學(xué)模型.運(yùn)用數(shù)值計(jì)算結(jié)果對(duì)拖纜系統(tǒng)動(dòng)力學(xué)仿真計(jì)算的穩(wěn)態(tài)構(gòu)型進(jìn)行了驗(yàn)證說(shuō)明，表明研究方法的正確性.

(2)拖纜的垂直長(zhǎng)度受飛機(jī)和拖體等參數(shù)的影響.飛機(jī)速度對(duì)系統(tǒng)工作的影響最大；為確保垂直度在要求范圍內(nèi)，在一定條件下，應(yīng)合理減小速度來(lái)增大拖纜垂直度；而增加末端拖體質(zhì)量類似于減少空氣阻尼的效果.

(3)滿足工程實(shí)際要求，需確保垂直度始終大于70%，系統(tǒng)正常工作的飛機(jī)速度范圍是v≤83 m·s-1.

[1] 李紅泉. 新型航空拖纜制作與檢測(cè)方法[J]. 兵工自動(dòng)化, 2014, 33(5): 19.

LI Hongquan. Manufacture and test method for a new towing cable[J]. Ordnance Industry Automation, 2014, 33(5): 19.

[2] 韓維，侯志強(qiáng)．機(jī)載拖曳系統(tǒng)的建模和若干動(dòng)力學(xué)問(wèn)題[C]∥中國(guó)2005年飛行力學(xué)與飛行實(shí)驗(yàn)學(xué)術(shù)年會(huì)論文集. 西安：飛行力學(xué)雜志社, 2005: 1-5．

HAN Wei, HOU Zhiqiang. The modeling and some dynamic problem of the trailing cable on an aircraft[C]∥The Flight Dynamic and Flight Experimentation Conference of China. Xi’an: The Flight Dynamic Press, 2005: 1-5.

[3] CAUGHEY T K. Whirling of a heavy chain[C]∥Proceedings of the 3rd U.S. National Congress of Applied Mechanics, American Society of Mechanical Engineers. New York: American Society of Mechanical Engineers, 1958: 101-108.

[4] LEMON G, FRASER W B. Steady-state bifurcations and dynamical stability of a heavy whirling cable acted on by aerodynamic drag[J]. Proceedings of the Royal Society of London A, 2001, 457(2009): 1021.

[5] SKOP R A, CHOO Y I. The Configuration of a cable towed in a circular path[J]. Journal of Aircraft, 1971, 8(11): 856.

[6] HUANG S L. Mathematical model for long cable by orbiting aircraft[R]. Johnsville: Naval Air Development Center, 1969.

[7] ZHU F, RAHN C D. Stability analysis of a circularly towed cable-body system[J]. Journal of Sound and Vibration, 1998, 217(3): 435.

[8] 鄭小洪，韓維，賈忠湖，等. 機(jī)載拖曳天線動(dòng)力學(xué)建模與仿真[J]. 飛行力學(xué), 2013, 31(3): 234.

ZHENG Xiaohong, HAN Wei, JIA Zhonghu,etal. Dynamic modeling and simulation of airborne trailing antenna[J]. Flight Dynamics, 2013, 31(3): 234.

[9] CHOO Y I, CASARELLA M J. A survey of analytical methods for dynamic simulation of cable-body systems[J]. Journal of Hydronautics, 1973, 7(4): 137.

[10] RUSSELL J J, ANDERSON W J. Equilibrium and stability of a whirling rod-mass system[J]. International Journal of Non-Linear Mechanics, 1977, 12(2): 91.

[11] JAMES M C, LOUIS V S, Thomas D S. Dynamic modeling of a trailing wire towed by an orbiting aircraft [J]. Journal of Guidance, Control, and Dynamics, 1995, 18(4): 875.

[12] PAUL W, PAVEL T. A study on the transitional dynamics of a towed-circular aerial cable system[C/CD] ∥ AIAA Atmospheric Flight Mechanics Conference and Exhibit. California: American Institute of Aeronautics and Astronautics, 2005.

[13] HOERNER S F. Fluid dynamic drag[M]. Bakersfield: Hoerner Fluid Dynamics, 1965.

[14] QUISENBERRY J, ARENA A. Dynamic simulation of low altitude aerial tow systems[C]∥AIAA Atmospheric Flight Mechanics Conference and Exhibit. Providence: AIAA, 2004: 4804-4813.

[15] BORST R G, GREISZ G F, QUYNN A G. Fuzzy logic control algorithm for suppressing E-6A long trailing wire antenna wind shear induced oscillations[R]. Washington: AIAA, 1993.