余本嵩，湯毓寧

（南京航空航天大學(xué)機(jī)械結(jié)構(gòu)力學(xué)及控制國(guó)家重點(diǎn)實(shí)驗(yàn)室，江蘇 南京 210016）

引 言

空間繩系系統(tǒng)以其低能耗、可重構(gòu)、易循環(huán)使用等優(yōu)勢(shì)已成為當(dāng)今航天領(lǐng)域探索的熱門課題之一，其在碎片捕獲、返回艙再入、廢棄物離軌、對(duì)地偵查、星際探測(cè)等諸多方面的不斷發(fā)展，更是使國(guó)內(nèi)外廣大學(xué)者在其理論、實(shí)驗(yàn)、應(yīng)用等方向投入大量精力［1-7］。

非軌道面繩系系統(tǒng)動(dòng)力學(xué)研究是其由理論到實(shí)際應(yīng)用的重要環(huán)節(jié)，已有相關(guān)科研成果呈現(xiàn)。比如，對(duì)于面內(nèi)外運(yùn)動(dòng)耦合的空間電動(dòng)力繩，Kojima 等［8］為了維持系統(tǒng)在傾斜軌道下的穩(wěn)態(tài)運(yùn)動(dòng)，設(shè)計(jì)了一套時(shí)滯反饋控制律，實(shí)現(xiàn)了對(duì)不規(guī)則振蕩的有效抑制。Ellis 等［9］理論分析自旋繩系衛(wèi)星時(shí)發(fā)現(xiàn)，微幅面外振動(dòng)不會(huì)對(duì)系統(tǒng)面內(nèi)運(yùn)動(dòng)產(chǎn)生影響，且相對(duì)高速的面內(nèi)自旋會(huì)使面外運(yùn)動(dòng)更加穩(wěn)定。Avanzini等［10］對(duì)非軌道面、雙金字塔型等繩系編隊(duì)進(jìn)行數(shù)值研究，算例表明系統(tǒng)穩(wěn)定性對(duì)系繩線密度、阻尼、剛度等參數(shù)較為敏感。楊育偉等［11］利用攝動(dòng)法解析得到面外電動(dòng)力繩系衛(wèi)星的近似周期解，并通過(guò)Poincaré 映射判定該周期運(yùn)動(dòng)的穩(wěn)定性。Paul 等［12］基于反饋線性化方法嘗試對(duì)繩系衛(wèi)星系統(tǒng)面內(nèi)/面外振蕩進(jìn)行控制，面內(nèi)振蕩得到有效抑制，面外振動(dòng)控制效果不理想。Yoon 等［13］通過(guò)等效縮比地面實(shí)驗(yàn)，研究空間系繩面外圓周運(yùn)動(dòng)，對(duì)系統(tǒng)振動(dòng)頻率等進(jìn)行討論，實(shí)驗(yàn)與有限元計(jì)算結(jié)果一致。在空間碎片移除應(yīng)用背景下，Cui 等［14］充分分析系繩纏繞等不穩(wěn)定因素，計(jì)算出初始面內(nèi)/面外傾角的合理范圍，以保證電動(dòng)力繩安全釋放。

對(duì)于包括混沌在內(nèi)的空間系繩非線性特性研究，仍主要聚焦于系統(tǒng)面內(nèi)運(yùn)動(dòng)。譬如，Nakanishi等［15］基于啞鈴模型研究軌道平面內(nèi)兩體繩系衛(wèi)星，數(shù)值結(jié)果表明，當(dāng)軌道偏心率大于0.3138 時(shí)系統(tǒng)面內(nèi)俯仰運(yùn)動(dòng)會(huì)出現(xiàn)混沌現(xiàn)象，并通過(guò)將不同軌道下面內(nèi)運(yùn)動(dòng)軌跡投影于van der Pol 平面以展示系統(tǒng)概周期、混沌特性。考慮主航天器剛體姿態(tài)，龐兆君等［16］開(kāi)展繩系航天器動(dòng)力學(xué)地面等效實(shí)驗(yàn)，依據(jù)天-地動(dòng)力學(xué)相似原理，揭示在軌剛體混沌運(yùn)動(dòng)，同時(shí)提出通過(guò)阻尼力矩對(duì)混沌進(jìn)行抑制。采用離散珠點(diǎn)模型，并計(jì)入熱沖擊、大氣阻尼等環(huán)境攝動(dòng)，Yu 等［17］討論了由軌道偏心率變化引起的繩系衛(wèi)星面內(nèi)俯仰運(yùn)動(dòng)分岔行為。Malashin 等［18］用非線性偏微分方程構(gòu)建空間繩系系統(tǒng)模型，并基于該模型提出邊界控制問(wèn)題，從而解析分析柔性系繩的橫向/縱向非線性振蕩及抑制方法。Lim 等［19］采用絕對(duì)節(jié)點(diǎn)坐標(biāo)法數(shù)值研究空間繩-網(wǎng)系統(tǒng)對(duì)軌道碎片的捕獲過(guò)程，提出大量關(guān)于柔性系繩纏繞、沖擊等非線性問(wèn)題。Aslanov等［20］分析一類面內(nèi)繩系拖曳系統(tǒng)，仿真結(jié)果顯示，隨著系繩長(zhǎng)度增加，面內(nèi)俯仰角及剛體姿態(tài)角的平衡位置皆會(huì)發(fā)生超臨界叉型分岔。Liu 等［21］探討在地-月引力場(chǎng)作用下含阻尼繩系系統(tǒng)的三體問(wèn)題，算例展示系統(tǒng)在L1和L2點(diǎn)附近存在概周期及混沌等運(yùn)動(dòng)。針對(duì)電動(dòng)力繩系航天器升軌技術(shù)，Li 等［22］通過(guò)模擬發(fā)現(xiàn)，無(wú)控狀態(tài)下系繩長(zhǎng)度、系統(tǒng)面內(nèi)外運(yùn)動(dòng)皆為非周期振蕩。Salazar 等［23］以數(shù)值方式呈現(xiàn)狀態(tài)保持階段面內(nèi)外耦合（即四維狀態(tài)空間）繩系衛(wèi)星的混沌運(yùn)動(dòng)，并基于Kelvin-Talt-Chetayev 定理設(shè)計(jì)拉力控制器對(duì)混沌進(jìn)行抑制。

可以看出，非軌道面繩系衛(wèi)星的周期運(yùn)動(dòng)穩(wěn)定性及控制方法得到了大量關(guān)注，但對(duì)其混沌動(dòng)力學(xué)卻鮮有涉及，而對(duì)面內(nèi)系統(tǒng)的混沌特性也多局限于數(shù)值討論，理論研究不充分。本文致力于面外繩系系統(tǒng)的混沌運(yùn)動(dòng)預(yù)測(cè)分析，具體討論面內(nèi)俯仰角保持恒定僅做面外滾轉(zhuǎn)的繩系衛(wèi)星，揭示面外攝動(dòng)引起的混沌，并給出混沌發(fā)生的必要條件，最后通過(guò)數(shù)值仿真進(jìn)行驗(yàn)證。

1 非軌道面滾轉(zhuǎn)繩系衛(wèi)星模型

由主星M、子星S 及起連接作用的絕熱系繩構(gòu)成的在軌繩系衛(wèi)星系統(tǒng)，如圖1所示。主星與子星皆為圓柱體，質(zhì)量分別為mM和mS?？紤]到狀態(tài)保持階段系繩始終處于緊繃狀態(tài)，故將其簡(jiǎn)化為線密度為ρt、長(zhǎng)度為l、直徑為dt的剛性桿，其兩端延長(zhǎng)線與衛(wèi)星圓柱體軸線重合。假設(shè)地球?yàn)榫|(zhì)球體，系統(tǒng)質(zhì)心o位于距地高度為Ho的繞地圓周軌道，在軌道平面Π 和非軌道平面Σ，系統(tǒng)運(yùn)動(dòng)可分別以面內(nèi)俯仰角θ和面外滾轉(zhuǎn)角?表示。為討論非軌道平面系統(tǒng)滾轉(zhuǎn)運(yùn)動(dòng)，施加與非軌道平面Σ 正交的控制力，使系統(tǒng)保持恒定的面內(nèi)俯仰角θC。此外，以ν表示繞地真近點(diǎn)角，以i表示軌道與赤道平面夾角。構(gòu)建固結(jié)于系統(tǒng)質(zhì)心o的軌道坐標(biāo)系o-xyz，其x軸指向系統(tǒng)運(yùn)動(dòng)方向，y軸正交于軌道平面，z軸指向地球質(zhì)心。

圖1 面外滾轉(zhuǎn)繩系衛(wèi)星系統(tǒng)Fig.1 Tethered satellite system with out-of-plane roll motion

基于第二類Lagrange 方程建立系統(tǒng)動(dòng)力學(xué)模型。首先，系統(tǒng)質(zhì)心繞地和其自身滾轉(zhuǎn)而產(chǎn)生的動(dòng)能可分別表示為：

式中 符號(hào)“′”表示變量對(duì)時(shí)間t求導(dǎo)，Ro表示系統(tǒng)質(zhì)心距地球質(zhì)心的位置標(biāo)量，=mM+mS+mt和為導(dǎo)出質(zhì)量參數(shù)，mt=ρtl為系繩質(zhì)量。同時(shí)，以無(wú)窮遠(yuǎn)處為勢(shì)能零點(diǎn)，列出系統(tǒng)勢(shì)能：

式中μE表示地球引力參數(shù)。

選取面外滾轉(zhuǎn)角?為廣義坐標(biāo)，將系統(tǒng)動(dòng)能T=Tt+Tr及以上勢(shì)能表達(dá)式代入第二類拉氏方程可得［24］：

式中Q?表示與廣義坐標(biāo)?對(duì)應(yīng)的廣義力。另外，圖1所示的真近點(diǎn)角ν，其相比于時(shí)間t能夠更直觀地展示系統(tǒng)軌道歷程，故以真近點(diǎn)角ν作為無(wú)量綱時(shí)間，用Lr表示參考長(zhǎng)度，引入無(wú)量綱變換：

則系統(tǒng)動(dòng)力學(xué)方程（3）可轉(zhuǎn)化為以下無(wú)量綱形式：

式中 符號(hào)“?”表示變量對(duì)無(wú)量綱時(shí)間ν求導(dǎo)。以上動(dòng)力學(xué)模型能夠?qū)哂泻愣鎯?nèi)俯仰角的繩系衛(wèi)星滾轉(zhuǎn)運(yùn)動(dòng)進(jìn)行描述。

2 攝動(dòng)激勵(lì)

考察大氣阻尼及太陽(yáng)光壓對(duì)系統(tǒng)面外滾轉(zhuǎn)的攝動(dòng)作用，基于虛功原理和力矩表達(dá)式可得到它們引起的攝動(dòng)力矩，即與面外滾轉(zhuǎn)角對(duì)應(yīng)的廣義力［25］：

式中 dF?p，t(s)為作用于系繩微元段ds的面外攝動(dòng)力，而由虛面外滾轉(zhuǎn)角δ?引起的微元段ds的虛位移表示為：

式中 ex，ey和ez分別表示與坐標(biāo)軸x，y和z對(duì)應(yīng)的單位矢量；s表示一個(gè)度量由子星指向主星的局部長(zhǎng)度坐標(biāo)，χ=(2mM+mt) (2mˉ)為導(dǎo)出參數(shù)。

作用于衛(wèi)星剛體和系繩微元段ds的大氣阻尼力為：

式中Cd，M(S)和Cd，t表示衛(wèi)星剛體和系繩的大氣阻尼系數(shù)，ρa(bǔ)，M(S)和ρa(bǔ)，t表示衛(wèi)星剛體和系繩所處軌道高度的大氣密度，V?r表示系統(tǒng)在面外方向相對(duì)于大氣的速度矢量，AM(S)表示衛(wèi)星剛體的迎風(fēng)面積。將衛(wèi)星剛體承受的大氣阻尼力對(duì)系統(tǒng)質(zhì)心取矩，而將式（8）中的第2 式代入式（6），得到大氣阻尼在廣義坐標(biāo)?方向作用于主衛(wèi)星剛體、子衛(wèi)星剛體及系繩的攝動(dòng)力矩：

式中ωE為地球自旋角速度。

系統(tǒng)處于始終向陽(yáng)的繞地軌道，則作用在衛(wèi)星剛體和系繩微元段ds的太陽(yáng)光壓力表示為：

式中Is=1372 W ?m?2表示太陽(yáng)在地球表面的入射輻射能流密度，vlc表示光速，β為太陽(yáng)光線與系統(tǒng)的夾角，Csr和Cdr表示系統(tǒng)表面的鏡面反射和漫反射系數(shù)，un和ut表示系統(tǒng)表面法向和光線投影方向的單位矢量。同樣方法，可得到太陽(yáng)光壓力在廣義坐標(biāo)?方向作用于主衛(wèi)星剛體、子衛(wèi)星剛體及系繩的攝動(dòng)力矩：

故相應(yīng)的由太陽(yáng)光壓引起的廣義力Q?sr=為：

將式（11）和（14）代入廣義力表達(dá)式Q?=Qd?+，則系統(tǒng)無(wú)量綱動(dòng)力學(xué)方程（5）可改寫為：

式中

至此，完成了空間環(huán)境攝動(dòng)激勵(lì)下非軌道面滾轉(zhuǎn)繩系衛(wèi)星系統(tǒng)模型的構(gòu)建。顯然，這是一個(gè)非自治的單自由度非線性系統(tǒng)。

3 混沌的必要條件

設(shè)向量?=[?1?2]T=[]T，則系統(tǒng)無(wú)量綱動(dòng)力學(xué)方程（15）可改寫為以下?tīng)顟B(tài)方程形式：

其中向量場(chǎng)和攝動(dòng)項(xiàng)的具體表達(dá)式為：

不難看出，攝動(dòng)項(xiàng)g(?，ν)的周期為2π，即有g(shù)(?，ν)=g(?，ν+2π)。值得注意的是，當(dāng)系統(tǒng)不考慮攝動(dòng)激勵(lì)，即g(?，ν)=0 時(shí)，其為一個(gè)Hamilton 系統(tǒng)，能夠計(jì)算出系統(tǒng)存在鞍點(diǎn)(?π 2，0)和異宿軌道：

將以上異宿軌線代入以下Melnikov 函數(shù)：

有

易見(jiàn)，在周期時(shí)間內(nèi)ν0∈[0，2π]，若系統(tǒng)滿足條件：

則Melnikov 函數(shù)存在簡(jiǎn)單零點(diǎn)。此必要條件意味著在充分小攝動(dòng)下，因系統(tǒng)異宿軌道（19）破缺導(dǎo)致穩(wěn)定流形與不穩(wěn)定流形在鞍點(diǎn)附近橫截相交，這可能誘發(fā)Smale 馬蹄意義下的混沌運(yùn)動(dòng)。

4 數(shù)值算例

取以下參數(shù)對(duì)繩系衛(wèi)星面外滾轉(zhuǎn)運(yùn)動(dòng)進(jìn)行數(shù)值仿真，以揭示系統(tǒng)混沌動(dòng)力學(xué)特性。主星和子星剛體質(zhì)量分別為mΜ=50 kg 和mS=50 kg，系繩線密度、長(zhǎng)度和直徑分別為ρt=5×10?3kg/m，l=1 km 和dt=0.5×10?3m。設(shè)衛(wèi)星剛體與系繩的大氣阻尼系數(shù)皆為Cd，M(S)(t)=2.2，主星和子星剛體的迎風(fēng)面積分別為AM=1.0 m2和AS=0.3 m2，系統(tǒng)表面的鏡面反射系數(shù)和漫反射系數(shù)分別為Csr=0.8 和Cdr=0.2。

基于解析的混沌判別式（22），對(duì)恒定面內(nèi)俯仰角、軌道傾角、軌道高度構(gòu)成的有限域內(nèi)參數(shù)(θC，i，Ho)依次取值計(jì)算，分析系統(tǒng)是否可能發(fā)生混沌，得到如圖2（a）所示混沌域，只要參數(shù)在此域內(nèi)，系統(tǒng)將可能發(fā)生混沌。選取截面Γ1={(θC，i，Ho)|θC=0 }和Γ2={(θC，i，Ho)|θC=3π 8 }，揭示恒定面內(nèi)俯仰角θC對(duì)混沌域的影響，如圖2（b）所示，可以看出面內(nèi)俯仰角越大混沌域越小。選取截面Γ3={(θC，i，Ho)|i=π 8 }和Γ4={(θC，i，Ho)|i=3π 8 }，圖2（c）呈現(xiàn)了軌道傾角i對(duì)混沌域的影響，混沌域隨軌道傾角增大而縮小。同樣，圖2（d）為軌道高度Ho對(duì)混沌域的作用，隨著質(zhì)心軌道升高混沌域?qū)⒆兇蟆?/p>

進(jìn)一步，以子星質(zhì)量、大氣阻尼系數(shù)、漫反射系數(shù)為例，研究系統(tǒng)結(jié)構(gòu)參數(shù)、大氣阻尼、太陽(yáng)光壓對(duì)混沌域的影響。基于初始設(shè)定參數(shù)，分別將子星質(zhì)量改設(shè)為mS=40 kg、大氣阻尼系數(shù)改設(shè)為Cd，M(S)(t)=1.0、漫反射系數(shù)改設(shè)為Cdr=0.6，能夠看出圖3（a），（b），（c）與圖2（a）所示混沌域存在差異?？梢园l(fā)現(xiàn)，混沌域?qū)ψ有琴|(zhì)量極為敏感，大氣阻尼、太陽(yáng)光壓對(duì)混沌域的影響也不可忽略，這些參數(shù)對(duì)繩系衛(wèi)星的動(dòng)力學(xué)設(shè)計(jì)至關(guān)重要。

圖2 參數(shù)(θC，i，Ho)的混沌域Fig.2 Chaotic domain governed by parameters(θC，i，Ho)

圖3 系統(tǒng)參數(shù)對(duì)混沌域的影響Fig.3 Effect of system parameters on chaotic domain

令系統(tǒng)的恒定面內(nèi)俯仰角、軌道傾角和軌道高度分別為θC=π 4，i=π 8 和Ho=350 km，進(jìn)行動(dòng)力學(xué)模擬?？捎?jì)算出當(dāng)前系統(tǒng)參數(shù)位于混沌域內(nèi)，即滿足混沌判別式（22），有|γ2|=1.1321×10?6<2.8554×10?6。繩系衛(wèi)星系統(tǒng)的非軌道面滾轉(zhuǎn)運(yùn)動(dòng)數(shù)值結(jié)果如圖4所示。圖4（a）表示系統(tǒng)面外滾轉(zhuǎn)角的時(shí)間歷程，可以清楚地看到滾轉(zhuǎn)角呈現(xiàn)非周期變化，故這必然是一個(gè)不規(guī)則運(yùn)動(dòng)。圖4（b）為滾轉(zhuǎn)運(yùn)動(dòng)的Poincaré 截面，能夠觀測(cè)到在不穩(wěn)定鞍點(diǎn)附近有大量的橫截異宿點(diǎn)存在，這是異宿軌道破缺所致，極易導(dǎo)致混沌發(fā)生。圖4（c）為系統(tǒng)的功率譜密度，其在(0，0.2 Hz)范圍內(nèi)存在密集的功率譜（這意味系統(tǒng)遍歷無(wú)窮多個(gè)不穩(wěn)定周期軌道），這是混沌運(yùn)動(dòng)的重要特征。圖4（d）為系統(tǒng)的最大Lyapunov 指數(shù)，其隨無(wú)量綱時(shí)間變化且始終大于0。綜合以上結(jié)果可以判斷，當(dāng)前系統(tǒng)的面外滾轉(zhuǎn)運(yùn)動(dòng)具有非常明顯的混沌動(dòng)力學(xué)特性。

圖4 具有混沌特性的面外滾轉(zhuǎn)運(yùn)動(dòng)Fig.4 Out-of-plane roll motion with chaotic characteristics

作為對(duì)比算例，將系統(tǒng)質(zhì)心的軌道高度改設(shè)為Ho=250 km，其他參數(shù)保持不變，則系統(tǒng)參數(shù)不再滿足混沌判別式（22），即|γ2|=9.1034×10?6>6.4413×10?6。此時(shí)，系統(tǒng)的滾轉(zhuǎn)動(dòng)力學(xué)行為如圖5所示。從圖5（a）可以發(fā)現(xiàn)，面外滾轉(zhuǎn)角變化看似具有一定周期性，而圖5（b）中Poincaré 截面上的閉軌線則表明此滾轉(zhuǎn)運(yùn)動(dòng)具有概周期特性。在圖5（c）和5（d）中，功率譜密度的幾個(gè)孤立尖峰（這意味系統(tǒng)存在固有頻率比為無(wú)理數(shù)的幾個(gè)不穩(wěn)定周期運(yùn)動(dòng)）和系統(tǒng)最大Lyapunov 指數(shù)最終趨于0，進(jìn)一步說(shuō)明系統(tǒng)做概周期運(yùn)動(dòng)。在此算例中，混沌現(xiàn)象沒(méi)有出現(xiàn)，這與混沌判別式（22）的判定結(jié)果一致。

圖5 具有概周期特性的面外滾轉(zhuǎn)運(yùn)動(dòng)Fig.5 Out-of-plane roll motion with quasi-periodic characteristics

5 結(jié) 論

本文通過(guò)正交方向抑制面內(nèi)俯仰角以保持繩系衛(wèi)星面內(nèi)俯仰角恒定，研究環(huán)境激勵(lì)下系統(tǒng)面外滾轉(zhuǎn)非線性特性?；贛elnikov 方法，判定穩(wěn)定流形與不穩(wěn)定流形在鞍點(diǎn)附近是否可能橫截相交，推導(dǎo)出系統(tǒng)面外滾轉(zhuǎn)的混沌判別式。獲得具有動(dòng)力學(xué)指導(dǎo)意義的參數(shù)混沌域，其對(duì)子星質(zhì)量等系統(tǒng)參數(shù)具有敏感性，并在域內(nèi)外分別發(fā)現(xiàn)混沌、概周期現(xiàn)象。