王 欽,何星星,文援蘭

(1.國(guó)防科學(xué)技術(shù)大學(xué) 航天與材料工程學(xué)院,湖南 長(zhǎng)沙 410073;2.上海衛(wèi)星工程研究所,上海 200240)

0 引言

當(dāng)前,航天領(lǐng)域中大型撓性航天器所受的關(guān)注越來(lái)越廣泛,撓性附件的復(fù)雜度也越來(lái)越高,在精確控制星體姿態(tài)的同時(shí)有效抑制柔性附件的振動(dòng),成為目前研究的重要方向[1]。撓性附件的運(yùn)動(dòng)需用無(wú)限自由度的分布參數(shù)描述,而撓性運(yùn)動(dòng)和剛體的姿態(tài)運(yùn)動(dòng)相互耦合,用常微分方程和偏微分方程共同描述的是非線性、非定常系統(tǒng)。特別是當(dāng)附件結(jié)構(gòu)非常復(fù)雜時(shí),求解更難。柔性航天器動(dòng)力學(xué)建模中,合理選擇航天器數(shù)學(xué)建模方法和原理可簡(jiǎn)化推導(dǎo)過(guò)程并有助于提高模型的工程實(shí)用性,常用方法有離散坐標(biāo)法、分布參數(shù)法、混合坐標(biāo)法和正則模態(tài)法等[2]。工程中常用的參數(shù)法或有限元法,都是將無(wú)限自由度的系統(tǒng)轉(zhuǎn)為有限自由度系統(tǒng),再行求解[3]。

相關(guān)領(lǐng)域的研究取得了大量成果[4、5]。特別是隨著計(jì)算機(jī)性能的提升,采用混合坐標(biāo)系方法建立由剛體姿態(tài)坐標(biāo)描述中心剛體的運(yùn)動(dòng)以及由離散模態(tài)坐標(biāo)描述撓性附件的運(yùn)動(dòng),截去高階模態(tài)對(duì)姿態(tài)控制系統(tǒng)的影響,即可轉(zhuǎn)化為有限自由度對(duì)系統(tǒng)的描述,這既能準(zhǔn)確描述系統(tǒng)運(yùn)動(dòng)及動(dòng)力學(xué)特性,又便于設(shè)計(jì)與分析控制系統(tǒng)。本文基于結(jié)構(gòu)動(dòng)力學(xué)及剛體運(yùn)動(dòng)學(xué)和動(dòng)力學(xué),建立多平臺(tái)聯(lián)合仿真系統(tǒng),對(duì)航天器結(jié)構(gòu)-姿態(tài)全系統(tǒng)進(jìn)行性能了分析及評(píng)估。

1 理論建模

1.1 撓性附件振動(dòng)特性

對(duì)航天器撓性附件結(jié)構(gòu)動(dòng)力學(xué)系統(tǒng),由有限元方法的基本格式,可得系統(tǒng)的求解方程為

對(duì)式(1)的大型特征值問(wèn)題,用有限元法求解時(shí),常見(jiàn)的解法有矩陣反迭代法、子空間迭代法、里茲向量直接疊加法和Lanczos向量的直接疊加法等,可解得方程的特征值ωi及特征向量φi。在此基礎(chǔ)上,將主坐標(biāo)轉(zhuǎn)化為模態(tài)坐標(biāo),即引入a(t)=Φx(t)=。此處：Φ為標(biāo)準(zhǔn)模態(tài);n為模態(tài)階次;x(t)為模態(tài)坐標(biāo),且x(t)=[x1… xn]T。由此可得模態(tài)坐標(biāo)的動(dòng)力學(xué)方程為

式中：Ω=diag[(ω1)2(ω2)2… (ωn)2];R(t)為模態(tài)坐標(biāo)下的廣義力。對(duì)應(yīng)的初始條件轉(zhuǎn)為：x0=ΦTMa0,。此處：a0為a的初值。

一般工程結(jié)構(gòu)可視為結(jié)構(gòu)型阻尼,根據(jù)Φ的正交性,有

式中;ξi為結(jié)構(gòu)阻尼比。因此,此時(shí)原動(dòng)力學(xué)方程可變?yōu)榛ハ嗖获詈系亩A常微分方程

式中：i=1,2,…,n。求解式(4)即可獲得每個(gè)振型的響應(yīng),將其疊加后,則每個(gè)節(jié)點(diǎn)的位移響應(yīng)為

1.2 撓性附件-航天器姿態(tài)耦合動(dòng)力學(xué)

對(duì)復(fù)雜的撓性附件結(jié)構(gòu),本文建立集中參數(shù)模型,用有限元法對(duì)系統(tǒng)進(jìn)行離散化,建立有限自由度的動(dòng)力學(xué)模型。

1.2.1 衛(wèi)星運(yùn)動(dòng)線動(dòng)量方程

設(shè)衛(wèi)星由1個(gè)本體結(jié)構(gòu)B和N個(gè)撓性附件組成,其中第i個(gè)撓性附件為Bi,如圖1所示。建立慣性坐標(biāo)系O1-x1y1z1。令衛(wèi)星未變形時(shí)的整星質(zhì)心位于點(diǎn)O,在點(diǎn)O處建立本體坐標(biāo)系O-xyz;第i個(gè)撓性附件與衛(wèi)星本體的連接點(diǎn)位于點(diǎn)Oi,在點(diǎn)Oi處建立撓性附件坐標(biāo)系Oi-xiyizi,點(diǎn)Oi在O-xyz系的矢量為lpi。對(duì)第i個(gè)附件的任一質(zhì)量元(點(diǎn)K)mik,在Oi-xiyizi系中的矢量位置為rik,衛(wèi)星變形后的位移為uik。則第i個(gè)附件的任一質(zhì)量元mik在O1-x1y1z1系中的線速度可表示為

式中：wik為O-xyz系中衛(wèi)星變形后點(diǎn)K的矢量坐標(biāo);ω為衛(wèi)星的絕對(duì)速度。

圖1 衛(wèi)星撓性附件結(jié)構(gòu)Fig.1 Flexible appendages structure of satellite

設(shè)撓性變形為小量,變形后的位移函數(shù)亦為一階小量,并作部分線性化處理,有wik=lpi+rik+uik,。

令VT=,則整星系統(tǒng)的線動(dòng)量

式中：mT為整星系統(tǒng)質(zhì)量。因點(diǎn)O為未變形前整星的質(zhì)心,則=0,故P=mTVT+。

根據(jù)Newton第二定律,衛(wèi)星運(yùn)動(dòng)的線動(dòng)量方程為

式中：F為作用于衛(wèi)星的外力。

在撓性附件的有限元計(jì)算中,令求解得到的撓性附件坐標(biāo)系的結(jié)構(gòu)振型為Φik=[φi1(rik) φi2(rik) … φin(rik)],k=1,2,…,n;固有頻率為ωi1,ωi2,…,ωin。根據(jù)式(5),uik=Φikqi(t)。此處：qi(t)為附件結(jié)構(gòu)離散的模態(tài)坐標(biāo),且qi(t)=[qi1(t) qi2(t) … qin(t)]。

設(shè)Oi-xiyizi系至O-xyz系的轉(zhuǎn)換陣為T(mén)SiB,則線動(dòng)量方程可改為

令BTrani為第i個(gè)撓性附件在O-xyz系中的平動(dòng)耦合系數(shù),且(此處BTrani為3×N維矩陣)。則

1.2.2 衛(wèi)星角動(dòng)量方程

設(shè)系統(tǒng)相對(duì)點(diǎn)O的絕對(duì)角動(dòng)量為L(zhǎng),系統(tǒng)相對(duì)點(diǎn)O的力矩為MSat,則根據(jù)角動(dòng)量定律,有

式中：IT為衛(wèi)星未變形時(shí)的慣量陣。

將有限元法所得的振型函數(shù)代入式(11),并統(tǒng)一轉(zhuǎn)為O-xyz系,即

1.2.3 撓性附件運(yùn)動(dòng)方程

式中：Ψ(VT,ω)為僅含VT,ω的項(xiàng)。

仍用FEA結(jié)果uik=Φikq(t),并轉(zhuǎn)為O-xyz系,再由振型函數(shù)的歸一化性質(zhì),即(Φik)T(TSiB)TTSiBΦik=IN(此處：IN為單元陣),以及B Rot i,(B Rot i)T,附件動(dòng)能可簡(jiǎn)化為

第i個(gè)撓性附件的應(yīng)變能為

式中：Kik為第k個(gè)節(jié)點(diǎn)的剛度陣;(Λi)2=diag[(ωi1)2(ωi2)2… (ωiN)2]。

式中：i=1,2,…,k。

1.3 解耦分析

聯(lián)立式(10)、(13)、(17),可得經(jīng)線性化處理的衛(wèi)星動(dòng)力學(xué)方程組。若推廣為更一般情況,即

a)考慮撓性附件的結(jié)構(gòu)阻尼ζ=[ζ1ζ2…ζn]T時(shí),式(17)應(yīng)為

b)衛(wèi)星帶角動(dòng)量為h的偏置動(dòng)量輪時(shí),式(13)應(yīng)為

式中：

其中：IB-C1(mB)-1C2為等效慣量陣。式(21)表明：角加速度項(xiàng)已與線加速度項(xiàng)解耦。

2 聯(lián)合仿真方法

衛(wèi)星在軌運(yùn)行過(guò)程中,大型撓性附件振動(dòng)會(huì)使姿態(tài)動(dòng)力學(xué)與結(jié)構(gòu)動(dòng)力學(xué)形成較緊密的耦合。大型撓性結(jié)構(gòu)如太陽(yáng)翼、天線、機(jī)械臂等會(huì)因衛(wèi)星調(diào)姿引起激烈振動(dòng),姿控系統(tǒng)的活動(dòng)部件如飛輪等在起動(dòng)、停止和反向等過(guò)程中會(huì)對(duì)結(jié)構(gòu)產(chǎn)生擾動(dòng),而這些振動(dòng)又會(huì)影響衛(wèi)星姿態(tài)。傳統(tǒng)衛(wèi)星設(shè)計(jì)方法中往往忽略姿態(tài)動(dòng)力學(xué)與結(jié)構(gòu)動(dòng)力學(xué)的耦合效應(yīng),單獨(dú)針對(duì)衛(wèi)星姿態(tài)動(dòng)力學(xué)與結(jié)構(gòu)動(dòng)力學(xué)進(jìn)行分析和設(shè)計(jì),故由此獲得的方案在實(shí)際應(yīng)用時(shí)的結(jié)果常與設(shè)計(jì)預(yù)期不符。

本文基于結(jié)構(gòu)動(dòng)力學(xué)平臺(tái)和姿態(tài)動(dòng)力學(xué)開(kāi)發(fā)平臺(tái),實(shí)現(xiàn)對(duì)全系統(tǒng)進(jìn)行聯(lián)合仿真研究,其總體方案如圖2所示。撓性附件的振動(dòng)與衛(wèi)星姿態(tài)系統(tǒng)的耦合主要通過(guò)撓性附件振動(dòng)特性、衛(wèi)星本體的姿態(tài)調(diào)整時(shí)產(chǎn)生的加速度和推力力矩形成耦合效應(yīng)。此外,姿態(tài)調(diào)整時(shí)坐標(biāo)系的變化會(huì)改變附件振型和耦合系數(shù),增大了求解難度。為此,本文用數(shù)值模擬方法實(shí)現(xiàn)姿態(tài)控制系統(tǒng)與撓性附件結(jié)構(gòu)動(dòng)力學(xué)的同步仿真,姿態(tài)控制系統(tǒng)將姿態(tài)角、加速度和推力力矩等數(shù)據(jù)送至附件振動(dòng)特性的計(jì)算;結(jié)構(gòu)動(dòng)力學(xué)分析將算得的附件結(jié)構(gòu)固有頻率、振型等數(shù)據(jù)送至姿態(tài)控制仿真系統(tǒng),循環(huán)直至整個(gè)姿態(tài)系統(tǒng)穩(wěn)定。

圖2 聯(lián)合仿真方案Fig.2 Co-simulation program

通過(guò)結(jié)構(gòu)動(dòng)力學(xué)平臺(tái)ANSYS提供的APDL標(biāo)量參數(shù)和數(shù)組參數(shù)可方便地實(shí)現(xiàn)有限元分析全過(guò)程的參數(shù)化,形成參數(shù)化的批處理文件,完成各種參數(shù)化分析,并可極大地提高分析效率?；诙鄤傮w系統(tǒng)的衛(wèi)星姿態(tài)動(dòng)力學(xué)建模技術(shù)已相對(duì)成熟。在考慮撓性附件時(shí),常采用混合坐標(biāo)法建立撓性附件的振動(dòng)方程和中心剛體的動(dòng)力學(xué)方程,聯(lián)合仿真程序流程如圖3所示。

圖3 聯(lián)合仿真程序流程Fig.3 Program process of co-simulation

3 仿真結(jié)果與分析

3.1 撓性附件振動(dòng)特性

以某雙自旋衛(wèi)星為例,以其中一星載天線作為衛(wèi)星附件結(jié)構(gòu),如圖4所示。該天線結(jié)構(gòu)主要包括與衛(wèi)星本體焊接的支撐板、半球柱網(wǎng)面結(jié)構(gòu)體、信號(hào)發(fā)射外伸端及螺釘?shù)雀綄贅?gòu)件。有限元建模時(shí),忽略螺帽、螺釘?shù)刃⌒蜆?gòu)件,對(duì)模型作部分簡(jiǎn)化。因支撐板與衛(wèi)星本體相連,在分析星載天線的振動(dòng)特性時(shí),可視之為固支結(jié)構(gòu)。

圖4 天線結(jié)構(gòu)Fig.4 Antenna structure

對(duì)天線結(jié)構(gòu)進(jìn)行模態(tài)分析,所得1～6階的固有頻率相應(yīng)為11.37,14.81,21.26,21.37,28.94,29.24 Hz,其中前4階模態(tài)振型如圖5所示。

圖5 天線結(jié)構(gòu)1～4階振型Fig.5 First 4 modes of antenna structure

3.2 附件振動(dòng)與姿態(tài)動(dòng)力學(xué)耦合算例

仍以某雙自旋衛(wèi)星為研究對(duì)象,建立其姿態(tài)動(dòng)力學(xué)模型,其中包括姿態(tài)控制模塊、球-管式章動(dòng)阻尼器模塊,以及噴氣推力器模塊等。設(shè)初始時(shí)刻,平臺(tái)相對(duì)衛(wèi)星質(zhì)心慣性系姿態(tài)[ψ θ φ]=[0° 30° 30°],平臺(tái)角速率在O-xyz系中角速率ω=[0 0 0]rad/s,輸入附件天線結(jié)構(gòu)的振動(dòng)特性,聯(lián)合仿真所得衛(wèi)星本體姿態(tài)角及角速率變化分別如圖6、7所示。在姿態(tài)控制系統(tǒng)的作用下,自旋軸方向單位矢量端點(diǎn)在質(zhì)心慣性系中投影的變化軌跡如圖8所示。

仿真結(jié)果表明：在噴氣姿態(tài)控制和被動(dòng)章動(dòng)阻尼控制下,衛(wèi)星平臺(tái)姿態(tài)角趨于目標(biāo)姿態(tài)角,角速率趨于零。星載天線結(jié)構(gòu)作高頻微幅振動(dòng),自由端最大振幅約10-6m量級(jí),其振動(dòng)響應(yīng)對(duì)衛(wèi)星本體的影響較小。另由圖8可知：在姿態(tài)機(jī)動(dòng)和附件結(jié)構(gòu)振動(dòng)相互耦合的影響下,衛(wèi)星姿態(tài)控制系統(tǒng)仍能保持穩(wěn)定。

圖6 衛(wèi)星姿態(tài)角變化過(guò)程Fig.6 Variation of satellite attitudeangle

圖7 衛(wèi)星角速率變化過(guò)程Fig.7 Variation of satelliteangular rate

圖8 自旋軸方向端點(diǎn)在質(zhì)心慣性系中投影Fig.8 Projection with direction of spin axis on barycentric inertia coordinates

4 結(jié)束語(yǔ)

本文以某雙自旋衛(wèi)星星載天線為撓性附件,建立了基于有限元方法和剛體動(dòng)力學(xué)結(jié)合的航天器姿態(tài)動(dòng)力學(xué)模型,并實(shí)現(xiàn)了航天器結(jié)構(gòu)-姿態(tài)的動(dòng)力聯(lián)合仿真。研究結(jié)果表明：通過(guò)此聯(lián)合仿真方法能實(shí)現(xiàn)航天器姿態(tài)控制全系統(tǒng)的性能分析和評(píng)估,在工程應(yīng)用中有一定的參考價(jià)值。

[1]MODI V.Attitude dynamics of satellites with flexible appendages：a brief review[J].Journal of Spacecraft and Rockets.1974(11)：743-751.

[2]曲廣吉,程道生.復(fù)合柔性結(jié)構(gòu)航天器動(dòng)力學(xué)建模研究[J].中國(guó)工程科學(xué),1999,1(2)：52-56.

[3]屠善澄,陳義慶.衛(wèi)星姿態(tài)動(dòng)力學(xué)與控制(I)[M].北京：中國(guó)宇航出版社,1999.

[4]CHOONGSEOK O H.Vibration control of flexible spacecraft under attitude maneuver using adaptive controller[C]//AIAA Guidance.Navigation and control Conference and Exhibit.Keystone,CO：AIAA,2006-6182：1-20.

[5]蔣建平,李東旭.帶撓性附件航天器剛?cè)狁詈蟿?dòng)力學(xué)[J].上海航天,2005,22(5)：24-27.

[6]王勖成.有限單元法[M].北京：清華大學(xué)出版社,2005.