穆瑞楠，王藝睿，譚述君,4，吳志剛,4，齊朝暉

(1. 大連理工大學(xué)工業(yè)裝備結(jié)構(gòu)分析國家重點實驗室，大連 116024；2. 中國科學(xué)院國家空間科學(xué)中心，北京 101408；3. 中國科學(xué)院大學(xué)，北京 101407；4. 大連理工大學(xué)航空航天學(xué)院，大連 116024)

0 引 言

空間太陽能電站(SSPS, Space Solar Power Station)是于1968年由美國科學(xué)家Peter Glaser提出的一種超大空間結(jié)構(gòu)概念[1]。因其可實現(xiàn)連續(xù)工作、能量利用率高等諸多優(yōu)點，近年來越來越受到各國和各研究機(jī)構(gòu)的關(guān)注。不僅相繼提出了數(shù)十種構(gòu)型設(shè)計方案，在能量傳輸、在軌組裝及動力學(xué)分析與控制等各個相關(guān)領(lǐng)域開展了大量工作[2]。由于尺寸和剛度的限制，這類超大空間結(jié)構(gòu)很難在地面開展實驗驗證工作。因此，借助于仿真手段準(zhǔn)確建模并分析其在軌的動力學(xué)特性至關(guān)重要。

絕大多數(shù)空間太陽能電站構(gòu)型尺寸巨大[3]，通常設(shè)計為模塊組裝的方式在軌構(gòu)建，包括桿件、薄膜、接口等基本組件，其建模難度高，模型包含的單元節(jié)點數(shù)目龐大，分析計算耗時較長，且不利于進(jìn)行大量的參數(shù)分析。特別是對于其在軌的動力學(xué)分析，需要同時考慮其軌道運(yùn)動和姿態(tài)運(yùn)動，使得模型具有較高的非線性，進(jìn)一步加大求解難度并降低求解效率。因此，對于空間太陽能電站的在軌動力學(xué)分析，基于其結(jié)構(gòu)特點，將其簡化為易于建模分析，并且具有一致動力學(xué)特性的簡單結(jié)構(gòu)十分必要。由于模塊組裝的在軌構(gòu)建方式，其結(jié)構(gòu)具有周期特性，因此采用基于能量等效原理的連續(xù)體等效方法最為合適。這種方法基于結(jié)構(gòu)的周期性，不會因基本單元數(shù)量的擴(kuò)充而導(dǎo)致計算量大幅增加，又能較好刻畫結(jié)構(gòu)的動力學(xué)特性。基于能量等效原理的連續(xù)體等效建模方法已有效應(yīng)用于國際空間站主支撐桁架結(jié)構(gòu)[4]以及大型環(huán)形可展天線的周邊桁架結(jié)構(gòu)[5]。

對于空間太陽能電站這種具有超大尺寸和超高柔性的空間結(jié)構(gòu)，其尺寸達(dá)到公里量級，其軌道運(yùn)動、姿態(tài)運(yùn)動和結(jié)構(gòu)振動之間的耦合作用不可忽略。一些學(xué)者對空間太陽能電站的耦合動力學(xué)開展了研究。Wie[6]基于Abacus電站構(gòu)型，分析了環(huán)境干擾力作用下的軌道運(yùn)動與姿態(tài)運(yùn)動，提出了軌道-姿態(tài)耦合控制方案。吳志剛[7]基于太陽塔式構(gòu)型，考慮地球扁率的引力場，在地球同步拉普拉斯軌道上研究了高階重力和力矩對空間太陽能電站姿軌運(yùn)動的影響，發(fā)現(xiàn)高階力對軌道運(yùn)動影響較大，而對姿態(tài)運(yùn)動影響較小可忽略。同時，超大結(jié)構(gòu)的結(jié)構(gòu)振動引起了更多學(xué)者的關(guān)注。Malla[8]、Ishimura[9]和穆瑞楠[10]建立了軸向振動空間啞鈴的軌道-姿態(tài)-振動耦合動力學(xué)模型，分別討論了初始條件、系統(tǒng)參數(shù)以及結(jié)構(gòu)尺寸對空間啞鈴動力學(xué)特性的影響。在穆瑞楠的研究工作中特別提到結(jié)構(gòu)振動引起了姿態(tài)運(yùn)動周期的改變，在較大初始姿態(tài)角下可能引起姿態(tài)運(yùn)動的翻滾現(xiàn)象。Silva[11-12]和Liu[13]以柔性梁為研究對象，前者給出了柔性梁的非線性擾動方程，并用擾動分析方法討論了系統(tǒng)可能出現(xiàn)的共振現(xiàn)象；而后者考慮動力剛化效應(yīng)，提出了參數(shù)激勵模型(PEM)，該模型可以更加精確地描述超大柔性空間結(jié)構(gòu)在持續(xù)轉(zhuǎn)動時的結(jié)構(gòu)振動。魏乙[14]將JAXA繩系構(gòu)型簡化為集中質(zhì)量-繩索-板模型，基于絕對節(jié)點坐標(biāo)法建立了耦合動力學(xué)方程，提出了結(jié)合辛Runge-Kutta方法的微分-代數(shù)方程求解方法，分析了系統(tǒng)參數(shù)對耦合特性的影響并驗證了求解方法的有效性。

本文基于能量等效的原則，將空間太陽能電站等效為柔性梁模型，并考慮重力梯度影響，建立姿態(tài)運(yùn)動-結(jié)構(gòu)振動耦合動力學(xué)模型；提出Runge-Kutta法和Newmark法相結(jié)合的耦合動力學(xué)方程的高效求解方法；最后對比不同參數(shù)下的仿真結(jié)果，并對超大柔性梁模型的在軌耦合特性進(jìn)行討論。

1 連續(xù)體等效建模

基于能量等效的連續(xù)體等效方法既可應(yīng)用于桁架結(jié)構(gòu)，也可應(yīng)用于框架結(jié)構(gòu)，只需要結(jié)構(gòu)具有周期性。多旋轉(zhuǎn)關(guān)節(jié)空間太陽能電站(MR-SPS)構(gòu)型的支撐部分為桁架梁組成的平面框架結(jié)構(gòu)，可分解為周期框架模塊和周期桁架模塊。首先針對周期桁架模塊，將桁架梁等效為連續(xù)梁模型；然后針對周期框架模塊，從而得到整體的等效連續(xù)梁模型?；谀芰康刃г淼倪B續(xù)體等效建模的理論方法如下。

周期模塊的各根構(gòu)件簡化為桿或梁模型。定義構(gòu)件局部坐標(biāo)系原點位于構(gòu)件中心，其x1軸方向沿構(gòu)件方向。每根構(gòu)件的變形為沿x1軸的拉伸變形以及繞y1軸和z1軸的彎曲變形，而其運(yùn)動為質(zhì)心沿x1軸、y1軸和z1軸的平動與繞質(zhì)心轉(zhuǎn)動。周期模塊總的變形能和動能可表示為

(1)

(2)

將周期模塊內(nèi)各點位移和應(yīng)變在周期模塊整體的幾何中心處Taylor展開，則可將周期模塊任一橫截面上各點的位移和應(yīng)變用周期模塊中心處的位移和轉(zhuǎn)角及應(yīng)變和曲率表示，即

u=Γuuc+Γεεc

(3)

其中uc表示周期模塊中心處的位移與轉(zhuǎn)角，而εc表示其應(yīng)變與曲率，Γu和Γε分別為對應(yīng)的展開系數(shù)矩陣。將式(3)代入式(1)和式(2)即可得到以中心點位移和應(yīng)變描述的周期模塊的變形能和動能。

另一方面，可建立與模塊相同長度的等效梁模型的變形能及動能。將等效梁模型在模塊長度內(nèi)的位移及應(yīng)變近似為均勻，則以中心點位移和應(yīng)變描述的變形能和動能表達(dá)式分別為

(4)

(5)

其中L為模塊長度。

分別對比周期模塊和等效梁模型的變形能和動能表達(dá)式，可得結(jié)構(gòu)的等效關(guān)系為

(6)

(7)

2 在軌柔性梁的動力學(xué)建模

基于哈密頓原理，考慮重力梯度影響，建立柔性梁的姿態(tài)運(yùn)動和結(jié)構(gòu)振動的耦合動力學(xué)模型。將地球視為理想球體，用Oe表示地心，柔性梁在軌道平面內(nèi)繞地球運(yùn)動，如圖1所示。

圖1 柔性梁結(jié)構(gòu)在軌運(yùn)行示意圖Fig.1 Schematic diagram of flexible beam structure in orbit

以柔性梁幾何中心的運(yùn)動作為其軌道運(yùn)動，軌道半徑R由地心指向柔性梁幾何中心，指向近地點的長半軸方向與軌道半徑矢量方向之間的夾角θ為軌道角，并用ωo表示軌道角速度。選取柔性梁幾何中心O為原點，柔性梁幾何中心處切線方向為x軸方向，建立柔性梁的固聯(lián)浮動坐標(biāo)系，則軌道半徑矢量方向與浮動坐標(biāo)系x軸方向的夾角φ為面內(nèi)姿態(tài)角，梁上各點偏移浮動坐標(biāo)系x軸的橫向位移為彎曲變形。坐標(biāo)量ρP由幾何中心O到點P變形前位置；橫向變形量wP由點P變形前位置指向變形后位置。則梁上點P相對于地心的位置rP為

(8)

(9)

其中mP表示點P的處的結(jié)構(gòu)質(zhì)量，即線密度。

對點P的動能式(9)在結(jié)構(gòu)長度上積分，得到整體柔性梁的動能為

(10)

對于在軌柔性梁結(jié)構(gòu)而言，其勢能包含重力勢能和變形能兩部分。這里考慮的柔性梁為純彎曲梁模型，其應(yīng)變能與教材中給出的歐拉梁應(yīng)變能一致[15]。因此柔性梁的總勢能為

(11)

其中μ為地球引力常數(shù)，EI為柔性梁的彎曲剛度，r為梁上某點到地心的距離，表示為

r2=(R2+x2+w2+2xRcosφ-2wRsinφ)

(12)

無非保守力做功的Hamilton原理為

(13)

將柔性梁的動能式(10)和勢能式(11)代入式(13)，得到柔性梁的在軌姿態(tài)-振動耦合動力學(xué)方程為

(14)

(15)

注意到耦合動力學(xué)方程中與重力梯度有關(guān)的項含有(r-3)項，該非線性項難以用有限元方法處理，因此需要將其近似展開。由式(12)得(r-3)的二項式展開為

(16)

其中參數(shù)ε=x/R,表示結(jié)構(gòu)尺寸與軌值半徑之比上式近似保留至參數(shù)ε的二次項以及w的線性項。將式代入式和式得到近似耦合動力學(xué)方程，其中含有參數(shù)ε的一階項和二階項分別對應(yīng)于重力梯度力/力矩的一階項和二階項。

由于柔性梁的在軌耦合動力學(xué)方程既包括偏微分方程，又含有在結(jié)構(gòu)長度上的積分項，因此采用有限元方法將結(jié)構(gòu)離散，引入只與時間t有關(guān)的節(jié)點位移變量和只與位置變量x有關(guān)的形函數(shù)，從而消去方程中關(guān)于位置變量x的微分和積分。梁上某點彎曲變形重新表述為

w(x,t)=N(x)ae(t)

(17)

其中N(x)表示單元內(nèi)的形函數(shù)，ae(t)表示單元內(nèi)的節(jié)點位移，這里梁單元采用經(jīng)典的二節(jié)點Hermite單元?；谟邢拊碚揫15]，利用梁橫向變形的新形式，可將在軌耦合動力學(xué)方程式(14)和式(13)分別轉(zhuǎn)化為

(18)

(19)

其中J為柔性梁面內(nèi)轉(zhuǎn)動慣量，M和K分別為柔性梁質(zhì)量陣和剛度陣，TG1和TG2分別為重力梯度力矩的一階項和二階項，PG1和PG2分別為重力梯度橫向分力的一階項和二階項。各項的表達(dá)式如下

其中Φ0、Φ1和Φ2是用于簡化標(biāo)記的矩陣和向量，其表達(dá)式為

3 模型參數(shù)及求解方法

多旋轉(zhuǎn)關(guān)節(jié)空間太陽能電站示意圖以及其中作支撐骨架的結(jié)構(gòu)力學(xué)模型如圖2(a)所示，其整體結(jié)構(gòu)由具有周期性的框架結(jié)構(gòu)組成，而周期框架結(jié)構(gòu)中的單根構(gòu)件又是由具有周期性的桁架結(jié)構(gòu)組成。桁架周期模塊和框架周期模塊以及它們的整體坐標(biāo)系如圖2(b)和(c)中藍(lán)色虛線(見電子版)框所示。

多旋轉(zhuǎn)關(guān)節(jié)空間太陽能電站構(gòu)型的幾何尺寸如表1所示，各桿件選用T700碳纖維復(fù)合材料，其材料參數(shù)在表1中給出。

利用基于能量等效的連續(xù)體等效方法，將該構(gòu)型等效為矩形等截面的均勻質(zhì)量柔性梁模型，其彎曲剛度EI= 1.56×109N·m2，密度與橫截面積的乘積(線密度)ρA= 6.20 kg/m。等效梁模型的結(jié)構(gòu)基頻為4.058 × 10-4Hz，而多旋轉(zhuǎn)關(guān)節(jié)構(gòu)型的有限元模型的基頻為3.841 × 10-4Hz，二者不僅具有相同量級，且差別僅為5%，能夠反應(yīng)其動力學(xué)現(xiàn)象。

圖2 桁架和框架周期模塊示意圖Fig.2 Diagram of truss and frame periodic modules

桁架模塊尺寸1 m×1 m×1 m框架模塊尺寸210 m×310 m總長度11.8 km桿件截面積2.359 cm2桿件彈性模量248 GPa桿件密度1750 kg/m3

涉及到的常數(shù)有：地球平均半徑RE=6378 km，地球引力常數(shù)為μ=3986 00 km3/s2。這里主要研究大尺寸結(jié)構(gòu)在低軌道(LEO)和地球同步軌道(GEO)下的動力學(xué)現(xiàn)象。選取柔性梁所處的低軌道的高度為322 km(即幾何中心軌道半徑R=6700 km)，此時的軌道頻率為1.833×10-4Hz，而地球同步軌道的高度為35 786 km(即幾何中心軌道半徑R=42 164 km)，此時的軌道頻率為1.161×10-5Hz. 同時，結(jié)構(gòu)尺寸是結(jié)構(gòu)基頻的重要影響因素，因此在50 m～10 000 m之間選取一系列值作為對比。初始姿態(tài)角對姿態(tài)運(yùn)動幅度有很大影響，這里取兩種情況：小姿態(tài)角(φ0=0.1 rad)和大姿態(tài)角(φ0=π/4 rad)；其他初始條件均設(shè)置為零。

在求解耦合方程式和式時，由于其已被轉(zhuǎn)化為二階常微分方程組，因此可采用經(jīng)典的Runge-Kutta 法求解。但由于有限元離散后的節(jié)點較多，使得方程維數(shù)較大，因此直接采用這種方法求解效率低。由于姿態(tài)和振動方程中系數(shù)相互包含，耦合求解時具有強(qiáng)非線性，因此難以使用常用的有限元方程求解方法Newmark法。對于非線性微分代數(shù)方程缺乏有效且高效的求解方法是超大柔性空間結(jié)構(gòu)耦合動力學(xué)研究所面臨的問題之一。本文將兩種方法結(jié)合，提出了Runge-Kutta 法和Newmark法同時迭代的方法，具體過程如圖3所示。與應(yīng)用Runge-Kutta法求解耦合方程相比，本文提出的方法在保證求解精度的同時效率大幅提升。在低軌道、小初始姿態(tài)角、結(jié)構(gòu)尺寸為1000 m的條件下，仿真時長1000 min，時間步長為0.1 s，本文提出的方法在效率上提高約752倍，兩種方法得出的響應(yīng)相差5.36%. 后續(xù)仿真結(jié)果均基于本文提出的改進(jìn)方法得出，其中時間步長取為0.1 s。

圖3 改進(jìn)的耦合模型求解算法流程圖Fig.3 Flow chart of improved algorithm for coupling model

在利用提出的方法求解動力學(xué)響應(yīng)時，姿態(tài)運(yùn)動方程式(18)中含有未知的初始加速度響應(yīng)，無法開始迭代。因此利用振動方程式(19)簡化姿態(tài)方程得

(20)

其中

4 姿態(tài)-振動耦合特性

4.1 結(jié)構(gòu)振動幅值影響分析

從結(jié)構(gòu)振動方程式(19)可以看出，姿態(tài)運(yùn)動對結(jié)構(gòu)振動既有慣性力的直接影響，也有改變重力梯度分力的間接影響。若忽略結(jié)構(gòu)變形對姿態(tài)運(yùn)動的影響，則姿態(tài)運(yùn)動方程式(18)可簡化為

(21)

將式(21)代入結(jié)構(gòu)振動方程(19)中發(fā)現(xiàn)，方程右端第一項(角加速度項)與第二項(重力梯度一階項)的一部分相互抵消，則結(jié)構(gòu)振動方程整理為

(22)

由上式可以看出，重力梯度的二階項PG2是激發(fā)結(jié)構(gòu)振動的主要因素，姿態(tài)運(yùn)動通過重力梯度的二階項PG2對結(jié)構(gòu)振動的間接影響。圖4分別給出了不同軌道下的大幅度姿態(tài)運(yùn)動下的結(jié)構(gòu)振動響應(yīng)。由重力梯度的二階項PG2的表達(dá)式可知，該項大小隨軌道半徑的增加反比例變化，同時以軌道角速度的平方量級減小，因此不同軌道下的振動量級近似相差3個量級。

圖4 大幅度姿態(tài)運(yùn)動下的結(jié)構(gòu)振動響應(yīng)Fig.4 Structural vibration response under attitude motion with large magnitude

由前文中的重力梯度的二階項PG2表達(dá)式可以看出，該項與質(zhì)量陣的比值以結(jié)構(gòu)尺寸的二次方量級增加。同時，剛度陣與質(zhì)量陣的比值以結(jié)構(gòu)尺寸的四次方量級減小。因此當(dāng)結(jié)構(gòu)基頻遠(yuǎn)大于姿態(tài)運(yùn)動頻率時，振幅以結(jié)構(gòu)尺寸的6次方量級增加。

圖5 結(jié)構(gòu)尺寸與振動幅度變化關(guān)系(LEO)Fig.5 Relationship between structural vibration magnitude and structure size(LEO)

圖5給出了低軌道下的結(jié)構(gòu)振幅隨結(jié)構(gòu)尺寸變化結(jié)果，在尺寸較小時，結(jié)構(gòu)振幅量級變化規(guī)律與理論結(jié)果一致，即按結(jié)構(gòu)尺寸的6次方量級增加。但仿真發(fā)現(xiàn)，結(jié)構(gòu)尺寸過大會引起系統(tǒng)不穩(wěn)定現(xiàn)象。以本文的參數(shù)為例，低軌道下結(jié)構(gòu)尺寸接近9450 m時，結(jié)構(gòu)振幅增長速度加快，最終系統(tǒng)運(yùn)動發(fā)散。這種不穩(wěn)定現(xiàn)象對應(yīng)的臨界尺寸與結(jié)構(gòu)柔性和軌道高度有關(guān)，軌道高度增加則對應(yīng)的臨界尺寸增大。

4.2 結(jié)構(gòu)振動頻率影響分析

由簡化的結(jié)構(gòu)振動方程式(22)可以得到結(jié)構(gòu)振動的第一階頻率為

(23)

其中ωs為結(jié)構(gòu)基頻。從式(23)可以看出，結(jié)構(gòu)振動頻率受到轉(zhuǎn)動耦合項(等式右端第二項)以及重力梯度項(等式右端第三項)的影響，使結(jié)構(gòu)振動周期出現(xiàn)波動現(xiàn)象。

考慮到在無控情況下姿態(tài)運(yùn)動接近于正弦振動，因此將姿態(tài)運(yùn)動近似表示為φ=φocos(κωot)，其中φ0為初始姿態(tài)角，kωo為姿態(tài)運(yùn)動頻率，k與φ0的取值相關(guān)，可通過與仿真得到的姿態(tài)運(yùn)動響應(yīng)擬合得到。將近似的姿態(tài)運(yùn)動代入式(23)并無量綱化得

1-3sin2(φ0coskωot)

(24)

注意到，式(24)中的轉(zhuǎn)動耦合項和重力梯度項中的軌道角速度處于三角函數(shù)內(nèi)，因此軌道角速度主要對頻率比的變化周期有影響，對幅值影響較小。

圖6 姿態(tài)運(yùn)動和重力梯度對振動頻率的影響(φ0 = π/4 rad)Fig.6 Effects of attitude motion and gravity gradient on frequency of structural vibration (φ0 = π/4 rad)

圖6給出了低軌道大幅姿態(tài)運(yùn)動下的轉(zhuǎn)動耦合項和重力梯度項對結(jié)構(gòu)振動頻率隨時間變化的影響曲線。對于大幅姿態(tài)運(yùn)動(φ0=π/4 rad)，對應(yīng)的k=1.5003，從圖6可以看出，重力梯度項波動幅度增大至一倍軌道頻率，時而使結(jié)構(gòu)振動頻率增大，時而使其減小。轉(zhuǎn)動耦合項仍使結(jié)構(gòu)振動頻率降低，其量級增大至兩倍軌道頻率。小幅姿態(tài)運(yùn)動下的結(jié)果類似，但幅值較小。

當(dāng)結(jié)構(gòu)基頻極低時，由式(23)可知結(jié)構(gòu)可能發(fā)生屈曲的失穩(wěn)現(xiàn)象，文獻(xiàn)[13]中對慣性穩(wěn)定的柔性梁的分析中也提及了這種現(xiàn)象。將令式(23)左側(cè)最小值等于0的結(jié)構(gòu)基頻稱為臨界結(jié)構(gòu)基頻。則對于不同幅度姿態(tài)運(yùn)動，對應(yīng)的臨界結(jié)構(gòu)基頻也不同。臨界基頻ωb滿足

1+3sin2(φ0coskωot)}

(25)

利用正弦函數(shù)的泰勒展開，以及二項式展開公式，式(25)右端括號內(nèi)的表達(dá)式可以轉(zhuǎn)化為參數(shù)λ的四階多項式形式。

(26)

其中λ=sinkωot。對式(20)求導(dǎo)并分析發(fā)現(xiàn)其導(dǎo)數(shù)表達(dá)式只有一個實根，同時該式的最高階的四次項系數(shù)為負(fù)，因此在全局有且只有一個極大值。對其導(dǎo)數(shù)表達(dá)式應(yīng)用一元三次方程求根公式，得到該極值點為

(27)

同時注意到λ在[-1,1]范圍內(nèi)取值，若極值點不在取值范圍內(nèi)則最大值點在邊界取值。經(jīng)數(shù)值驗證，在初始姿態(tài)角φ0在[0 , 0.476π]時，最大值在λ = -1處取得，則臨界振動頻率ωb滿足

(28)

圖7 臨界振動頻率隨初始姿態(tài)角變化曲線Fig.7 Variation of critical frequency of structural vibration with respect to initial attitude angle

圖7根據(jù)式(28)給出了低軌道和地球同步軌道下的臨界振動頻率隨初始姿態(tài)角變化曲線，發(fā)現(xiàn)不論在低軌道還是地球同步軌道，耦合作用對結(jié)構(gòu)振動頻率的影響都存在。初始姿態(tài)角在0.386π rad時臨界振動頻率最大，此時影響幅值約為2.11倍軌道頻率。因此，在設(shè)計梁式空間結(jié)構(gòu)時，考慮耦合作用對結(jié)構(gòu)振動頻率的影響，才能保證結(jié)構(gòu)不會發(fā)生屈曲不穩(wěn)定現(xiàn)象。

4.3 結(jié)構(gòu)振動對姿態(tài)運(yùn)動的影響

若引入小角度假設(shè)條件，并忽略結(jié)構(gòu)振動對轉(zhuǎn)動慣量的影響，則姿態(tài)運(yùn)動方程簡化為

(29)

其中

僅考慮柔性梁一階振型的振動，令

a=w0Ω(x)sin(ωst)

(30)

其中w0為最大振動幅值，Ω(x)為一階振型。將式(30)代入式(29)，得到姿態(tài)運(yùn)動有系統(tǒng)阻尼的自振頻率ωa為

(31)

其中阻尼比為

(32)

阻尼比的最大值可表示為

(33)

其中α是系統(tǒng)常數(shù)，與材料屬性、梁振型等參數(shù)有關(guān)。由式可以看出，姿態(tài)運(yùn)動頻率受結(jié)構(gòu)振動幅度和結(jié)構(gòu)尺寸的影響。結(jié)構(gòu)振動幅度增大時，姿態(tài)運(yùn)動頻率降低。同時注意到結(jié)構(gòu)振動幅度與結(jié)構(gòu)尺寸的6次方量級關(guān)系，因此結(jié)構(gòu)尺寸增加也導(dǎo)致姿態(tài)運(yùn)動頻率降低。

圖8分別給出了小幅度姿態(tài)運(yùn)動(φ0=0.1 rad)和大幅度姿態(tài)運(yùn)動(φ0= π/4 rad)下的剛體與柔性體姿態(tài)角響應(yīng)對比。為了顯示出柔性引起的姿態(tài)周期變化，圖中截取了仿真時間3200～3300 min的結(jié)果。從圖中可以發(fā)現(xiàn)，隨著仿真時間增加，柔性體姿態(tài)運(yùn)動越來越偏離剛體的姿態(tài)運(yùn)動，即柔性體姿態(tài)運(yùn)動周期較大。對比兩圖發(fā)現(xiàn)，初始姿態(tài)角越大，這種結(jié)構(gòu)柔性引起的姿態(tài)運(yùn)動周期增大現(xiàn)象越明顯，這是由于初始姿態(tài)角越大，姿態(tài)運(yùn)動引起的結(jié)構(gòu)振動幅度越大，反過來使得姿態(tài)運(yùn)動阻尼比越大，因此姿態(tài)運(yùn)動周期增大就越明顯。

圖8 柔性體與剛體的姿態(tài)角響應(yīng)比較(LEO)Fig.8 Comparison of attitude angle between flexible and rigid model(LEO)

圖9給出了地球同步軌道下的大幅度姿態(tài)運(yùn)動時的剛?cè)峤Y(jié)構(gòu)姿態(tài)角響應(yīng)結(jié)果對比，結(jié)構(gòu)柔性引起的姿態(tài)運(yùn)動周期增大現(xiàn)象仍然存在。分析式(33)可知，軌道角速度降低原本是會增大姿態(tài)運(yùn)動阻尼比，但是由圖4可知此時結(jié)構(gòu)振動幅度大幅下降，因此在地球同步軌道下剛體和柔性體的姿態(tài)運(yùn)動周期差別較小。

圖9 柔性體與剛體的姿態(tài)角響應(yīng)比較(GEO)Fig.9 Comparison of attitude angle between the flexible and rigid model(GEO)

5 結(jié)論與展望

本文將MR-SPS電站構(gòu)型等效為柔性梁模型，并考慮重力梯度影響，建立了姿態(tài)運(yùn)動和結(jié)構(gòu)振動耦合動力學(xué)模型，利用改進(jìn)算法分析了不同參數(shù)取值下兩者耦合關(guān)系。得到主要結(jié)論如下：1)所提出的基于Runge-Kutta 法和Newmark法的改進(jìn)算法適用于姿態(tài)運(yùn)動-結(jié)構(gòu)振動耦合方程求解，大幅提高了計算效率；2)在軌柔性梁模型中的重力梯度項至少近似保留至ε的二階展開項，尺寸過大可能導(dǎo)致結(jié)構(gòu)屈曲不穩(wěn)定；3)姿態(tài)運(yùn)動和結(jié)構(gòu)振動耦合效應(yīng)導(dǎo)致結(jié)構(gòu)振動頻率降低，姿態(tài)運(yùn)動周期增大，在GEO軌道仍有這種現(xiàn)象。為保證系統(tǒng)在軌穩(wěn)定性，超大型空間結(jié)構(gòu)設(shè)計時應(yīng)考慮耦合因素影響，并需要進(jìn)一步開展與結(jié)構(gòu)振動抑制協(xié)同的姿態(tài)控制策略研究。

參 考 文 獻(xiàn)

[1] Glaser P E. Power from the sun: its future[J]. Science, 1968, 162(3856):857-861.

[2] 侯欣賓, 王立. 空間太陽能電站技術(shù)發(fā)展現(xiàn)狀及展望[J]. 中國航天, 2015, (02):12-15. [Hou Xin-bin, Wang Li. Present situation and prospect for the solar power station [J]. Space System and Technology, 2015, (02):12-15.]

[3] 侯欣賓, 王立, 張興華等. 多旋轉(zhuǎn)關(guān)節(jié)空間太陽能電站概念方案設(shè)計[J]. 宇航學(xué)報, 2015, (11): 1332-1338. [Hou Xin-bin, Wang Li, Zhang Xing-hua, et al. Concept design on multi-rotary joints SPS [J]. Journal of Astronautics, 2015, (11): 1332-1338.]

[4] Noor A K, Mikulas M M. Continuum modeling of large lattice structures: status and projections[M]// Atluri S N, Anthony K A. Large space structures: dynamics and control. Heidelberg, Berlin: Springer, 1988:1-34.

[5] 劉福壽, 金棟平. 環(huán)形桁架結(jié)構(gòu)徑向振動的等效圓環(huán)模型[J]. 力學(xué)學(xué)報, 2016, 48(5):1184-1191. [Liu Fu-shou, Jin Dong-ping. Equivalent circular ring model for the radial vibration analysis of hoop truss structures[J]. Chinese Journal of Theoretical and Applied Mechanics, 2016, 48(5):1184 - 1191.]

[6] Wie B, Roithmayr C M. Attitude and orbit control of a very large geostationary solar power satellite[J]. Journal of Guidance, Control, and Dynamics. 2005, 28(3):439-451.

[7] 吳志剛, 劉玉亮, 張開明, 等. 高階重力和力矩對空間太陽能電站運(yùn)動的影響[J]. 空間控制技術(shù)與應(yīng)用, 2016, 42(4):1-5. [Wu Zhi-gang, Liu Yu-liang, Zhang Kai-ming, et al. The influences of higher order gravitational force and torque to the motion of space solar power station [J]. Aerospace Control and Application, 2016, 42(4):1-5.]

[8] Malla R B. Structural and orbital conditions on response of large space structures[J]. Journal of Aerospace Engineering, 1993, 6(2):115-132.

[9] Ishimura K, Higuchi K. Coupling among pitch motion, axial vibration, and orbital motion of large space structures[J]. Journal of Aerospace Engineering, 2008, 21(2):61-71.

[10] 穆瑞楠, 譚述君, 吳志剛. 重力梯度對超大柔性空間結(jié)構(gòu)在軌動力學(xué)特性的影響[J]. 國防科技大學(xué)學(xué)報, 2017, 39(3):7-14. [Mu Rui-nan, Tan Shu-jun, Wu Zhi-gang. Effect of gravity gradient on dynamical characteristics of very large flexible space structures in orbit[J]. Journal of National University of Defense Technology, 2017, 39(3):7- 14.]

[11] Dasilva M M, Zaretzky C L. Nonlinear dynamics of a flexible beam in a central gravitational field—I. equations of motion[J]. International Journal of Solids and Structures. 1993, 30(17):2287-2299.

[12] Dasilva M M, Zaretzky C L. Nonlinear dynamics of a flexible beam in a central gravitational field—II. nonlinear motions in circular orbit[J]. International Journal of Solids and Structures. 1993, 30(17): 2301-2316.

[13] Liu Y, Wu S, Zhang K, et al. Parametrical excitation model for rigid-flexible coupling system of solar power satellite[J]. Journal of Guidance, Control, and Dynamics. 2017:1-8.

[14] 魏乙, 鄧子辰, 李慶軍,等. 繩系空間太陽能電站動力學(xué)響應(yīng)分析[J]. 宇航學(xué)報, 2016, 37(9): 1041-1048. [Wei Yi, Deng Zi-Chen, Li Qing-jun, et al. Analysis of dynamic response of tethered space solar power station [J]. Journal of Astronautics, 2016, 37(9):1041-1048.]

[15] 王勖成. 有限單元法[M]. 北京: 清華大學(xué)出版社, 2003:306-311.