基于ANCF的薄膜太陽帆自旋展開動(dòng)力學(xué)模擬

薛鵬聰，劉 鋮,2,3，水小平

（1. 北京理工大學(xué) 宇航學(xué)院，北京 100081；2. 中物院高性能數(shù)值模擬軟件中心，北京 100088；3. 北京應(yīng)用物理與計(jì)算數(shù)學(xué)研究所，北京 100088）

引 言

膜結(jié)構(gòu)能夠在保證結(jié)構(gòu)可靠性的同時(shí)耗費(fèi)更少的材料包圍更大的空間，并具有自重更輕、造型優(yōu)美、便于拆卸等優(yōu)點(diǎn)，因此在航天工程領(lǐng)域得到了廣泛的關(guān)注，被期望作為執(zhí)行深空探測(cè)任務(wù)的新一代創(chuàng)新空間結(jié)構(gòu)，如薄膜太陽帆推進(jìn)系統(tǒng)、通訊衛(wèi)星中大型薄膜天線反射器及太陽能電站中電池陣列結(jié)構(gòu)等[1]。這類空間結(jié)構(gòu)在發(fā)射前通常處于折疊狀態(tài)，入軌后在驅(qū)動(dòng)裝置的作用下展開，并進(jìn)入工作狀態(tài)。然而對(duì)于在太空中展開太陽帆這類具有高收納比的大型薄膜結(jié)構(gòu)，無疑是一項(xiàng)復(fù)雜的操作，且展開后的有效面積及膜面精度均會(huì)對(duì)其所承擔(dān)的功能產(chǎn)生一定影響。因此，進(jìn)入太空中的薄膜航天器能否順利展開至理想工作狀態(tài)是這類航天器折展系統(tǒng)設(shè)計(jì)中的關(guān)鍵技術(shù)?？紤]到地面實(shí)驗(yàn)的成本過大，且難以針對(duì)大型薄膜結(jié)構(gòu)提供在軌微重力真空環(huán)境[2]。于是，采用高效的數(shù)值計(jì)算方法對(duì)空間帆膜結(jié)構(gòu)的動(dòng)力學(xué)展開過程進(jìn)行仿真模擬是此類薄膜航天器結(jié)構(gòu)設(shè)計(jì)的重要手段。Shirasawa等[3]采用多質(zhì)點(diǎn)模型對(duì)“伊卡洛斯”（Icarus）太陽帆展開動(dòng)力學(xué)進(jìn)行數(shù)值模擬，采用由彈簧和阻尼器連接的質(zhì)點(diǎn)代替?zhèn)鹘y(tǒng)薄膜單元，近似簡(jiǎn)化了膜結(jié)構(gòu)建模，降低了計(jì)算成本。Miyazaki[4]提出了一種剛度退化的褶皺模型，通過引入剛度折減系數(shù)描述薄膜壓縮時(shí)的本構(gòu)關(guān)系，建立了薄膜太陽帆的柔性多體系統(tǒng)有限元模型，并將其應(yīng)用于薄膜太陽帆的展開動(dòng)力學(xué)分析。

針對(duì)類似薄膜結(jié)構(gòu)的柔性多體系統(tǒng)的大轉(zhuǎn)動(dòng)、大變形問題，柔性多體系統(tǒng)動(dòng)力學(xué)領(lǐng)域提出了絕對(duì)節(jié)點(diǎn)坐標(biāo)方法（Absolute Nodal Coordinate Formulation，ANCF）[5]。基于ANCF方法，國內(nèi)學(xué)者針對(duì)薄膜太陽帆等多柔體航天器的展開動(dòng)力學(xué)進(jìn)行了較為深入的研究[6-9]。趙將等[7]采用ANCF薄板單元，對(duì)考慮粘彈性阻尼的“伊卡洛斯”薄膜太陽帆進(jìn)行了展開動(dòng)力學(xué)分析。Liu等[8]系統(tǒng)地對(duì)ANCF板/殼單元進(jìn)行了研究，結(jié)合SRM模型提出了縮減的ANCF薄膜單元，并通過薄膜太陽帆的自旋展開動(dòng)力學(xué)模擬驗(yàn)證了所提出方法的有效性。然而，由于薄膜結(jié)構(gòu)在宏觀響應(yīng)中拉壓模量呈現(xiàn)的巨大差異，系統(tǒng)非線性特性顯著，導(dǎo)致薄膜結(jié)構(gòu)在靜力學(xué)/動(dòng)力學(xué)模擬中極易發(fā)生計(jì)算不穩(wěn)定現(xiàn)象。

針對(duì)拉壓不同模量彈性結(jié)構(gòu)，前蘇聯(lián)學(xué)者阿姆巴book=218,ebook=122爾楚米揚(yáng)（Ambartsumyan）[10]系統(tǒng)地提出了不同模量彈性理論，并給出了其本構(gòu)關(guān)系的基本假定，推導(dǎo)了二維平面問題和三維空間問題的廣義彈性定律。之后許多研究者基于不同角度對(duì)Ambartsumyan的不同模量彈性理論進(jìn)行了發(fā)展和修正[11]。近年來，針對(duì)雙模量結(jié)構(gòu)的數(shù)值算法得到了發(fā)展，張亮等[12]基于參數(shù)變分原理和有限元方法，在主應(yīng)力方向上建立了統(tǒng)一的含參數(shù)變量的本構(gòu)方程，將不同模量問題演化為一個(gè)易于求解的標(biāo)準(zhǔn)互補(bǔ)二次規(guī)劃問題，避免了在有限元迭代過程中不斷根據(jù)應(yīng)力狀態(tài)更新剛度矩陣，并將其應(yīng)用于二維薄膜的起褶分析，呈現(xiàn)出不錯(cuò)的收斂效率。杜宗亮等[13]基于不同模量彈性理論推導(dǎo)了一種具有漸進(jìn)二次收斂效率的切線剛度矩陣，通過將其應(yīng)用于二維薄膜的應(yīng)力水平和起褶區(qū)域的預(yù)測(cè)，驗(yàn)證了該方法在處理不同模量薄膜結(jié)構(gòu)時(shí)的高效性、穩(wěn)定性。

本文旨在基于ANCF建模方法，整合張力場(chǎng)與不同模量彈性理論[10]，提出一種考慮不同模量的ANCF薄膜單元，以實(shí)現(xiàn)空間薄膜結(jié)構(gòu)的高效動(dòng)力學(xué)模擬。為提高單元在處理拉壓不同模量薄膜問題時(shí)的算法穩(wěn)定性，推導(dǎo)了單元精確的切線剛度矩陣，并通過若干靜力學(xué)算例驗(yàn)證了單元的高效性。此外，針對(duì)一正六邊形薄膜太陽帆，建立其展開過程多柔體系統(tǒng)動(dòng)力學(xué)模型，分析了薄膜太陽帆的結(jié)構(gòu)參數(shù)對(duì)其展開過程的影響，數(shù)值仿真結(jié)果為這類空間帆膜結(jié)構(gòu)的設(shè)計(jì)和優(yōu)化提供了理論指導(dǎo)。

1 基于ANCF描述的薄膜單元

如圖1所示，點(diǎn)P(ξ,η)為初始構(gòu)形中薄膜單元上一點(diǎn)，該點(diǎn)的位置矢量r0(ξ,η)表示為

book=219,ebook=123

2 基于不同模量彈性理論的彈性力及切線剛度

2.1 薄膜張力場(chǎng)理論與不同模量彈性理論

張力場(chǎng)理論[16]認(rèn)為薄膜的彎曲剛度為零，當(dāng)薄膜內(nèi)出現(xiàn)壓應(yīng)力時(shí)，一般通過面外變形即薄膜起褶的形式釋放壓應(yīng)力，張力場(chǎng)理論將其處理為薄膜的面內(nèi)收縮。起褶后的薄膜處于單軸拉伸應(yīng)力狀態(tài)，褶皺的方向平行于大主應(yīng)力的方向，且垂直于褶皺方向的小主應(yīng)力為零。

book=220,ebook=124

2.2 考慮不同模量的ANCF薄膜單元的切線剛度

3 基于廣義-α算法的多體系統(tǒng)動(dòng)力學(xué)方程

圖2給出了薄膜結(jié)構(gòu)的動(dòng)力學(xué)分析流程圖，在預(yù)處理模塊中，根據(jù)薄膜結(jié)構(gòu)的初始構(gòu)形以及材料參數(shù)和單元信息等建立了精確的ANCF薄膜單元，給出了系統(tǒng)的常數(shù)質(zhì)量矩陣。然后，采用高斯消元法消除冗余book=221,ebook=125約束，得到系統(tǒng)的稀疏矩陣。通過假定初始時(shí)刻系統(tǒng)的應(yīng)力場(chǎng)確定初始彈性矩陣，并組裝系統(tǒng)的彈性力和雅克比矩陣，然后代入到式（36）～（39）的廣義-α算法的求解過程中，根據(jù)求解結(jié)果更新系統(tǒng)位移場(chǎng)和彈性矩陣，如此循環(huán)迭代直至達(dá)到給定的收斂準(zhǔn)則，求解完成。

圖2 膜結(jié)構(gòu)動(dòng)力學(xué)分析計(jì)算方法Fig. 2 Dynamic analysis and calculation method of membrane structure

4 數(shù)值仿真驗(yàn)證

為便于區(qū)分，將式（29）對(duì)應(yīng)的線性近似的切線剛度模型記為ANCF LM，將式（30）對(duì)應(yīng)的切線剛度模型記為ANCF TM?；诒疚奶岢龅挠?jì)算方法，本節(jié)主要針對(duì)薄膜結(jié)構(gòu)的一些經(jīng)典案例進(jìn)行了考察。4.1節(jié)、4.2節(jié)和4.3節(jié)分別驗(yàn)證文所提出的計(jì)算框架在求解不同模量小變形、大變形，以及面外變形問題時(shí)的高效性和精確性。

4.1 基于ANCF的薄膜單元小變形分析

圖3所示為一個(gè)驗(yàn)證薄膜起褶模型各種數(shù)值算法的經(jīng)典算例，矩形薄膜寬度H，厚度t，其上下邊界施加有大小σ0的均勻預(yù)應(yīng)力，左右兩端施加有軸向集中力P=σ0tH，同時(shí)左右邊界均施加有集中力矩M。隨著力矩的增加，膜的下邊緣開始出現(xiàn)一定寬度的褶皺區(qū)域。考慮采用ANCF DM薄膜單元計(jì)算框架，對(duì)上述薄膜結(jié)構(gòu)進(jìn)行了數(shù)值模擬，仿真過程中，不考慮薄膜的壓縮模量。

圖3 純彎曲矩形薄膜受力示意圖Fig. 3 Stress diagram of Pure Bending Rectangular membrane

文獻(xiàn)[19]針對(duì)上述模型，采用解析方法對(duì)其起褶區(qū)域以及面內(nèi)變形進(jìn)行了研究，并得出一系列關(guān)系式。其中，式（40）給出了矩形膜起褶區(qū)域的寬度b與外力矩M的解析關(guān)系式

式（41）給出了不同工況下，薄膜右邊界上一點(diǎn)的柯西應(yīng)力與外力矩的解析關(guān)系式

式（42）給出了薄膜變形后曲率κ與外力矩M之間的解析關(guān)系式

圖4給出了起褶區(qū)域?qū)挾鹊慕馕鼋夂蛿?shù)值試驗(yàn)結(jié)果，不難看出即使在起褶區(qū)域超過90%時(shí)，數(shù)值模擬結(jié)果仍然具有很高的精度。

圖4 褶皺帶寬與外力矩的關(guān)系Fig. 4 Relationship between fold bandwidth and external torque

book=222,ebook=126圖5則給出了不同外力矩作用下，右邊界上一點(diǎn)柯西應(yīng)力的解析解與數(shù)值試驗(yàn)結(jié)果，兩者的吻合性比較好。

圖5 右邊界水平應(yīng)力與外力矩的關(guān)系Fig. 5 Relationship between horizontal stress at right boundary and external moment

圖6給出了薄膜結(jié)構(gòu)曲率的解析解和數(shù)值解隨外力矩的變化趨勢(shì)，分析得知，兩者的趨勢(shì)是非常一致的。

圖6 純彎曲矩形薄膜的彎矩–曲率關(guān)系Fig. 6 Moment-curvature relationship of rectangular membrane

通過上述各項(xiàng)數(shù)值試驗(yàn)結(jié)果與解析結(jié)果的對(duì)比可以得出，本文所提出的方法在研究薄膜起褶的問題時(shí)具有很高的精度。

4.2 基于ANCF的不同模量結(jié)構(gòu)大變形分析

如圖7所示為一L形梁結(jié)構(gòu)，其左端為固定約束，自由端受水平方向集中力F= 40 000 N。固定其受拉方向的彈性模量為Et= 3×107Pa，泊松比為vt= 0.3，依次分別取壓縮模量為Ec=3×107Pa、Ec=1.5×107Pa、Ec= 6 ×106Pa、Ec=3×106Pa 以及Ec= 1.5 × 106Pa 時(shí)的工況，對(duì)應(yīng)的拉壓比分別為Et/Ec=1、2、5、10、20。

為驗(yàn)證本文所提出計(jì)算方法的有效性，本節(jié)基于ABAQUS平臺(tái)的擴(kuò)展程序接口進(jìn)行了二次開發(fā)，并將UMAT的結(jié)果作為ANCF描述的不同模量問題的參考解。需要指出的是，UMAT這里采用CPS4I平面應(yīng)力四邊形單元和自適應(yīng)時(shí)間增量步進(jìn)行計(jì)算，在拉壓比較大的情況下需要較多的加載步。而本文所提出的ANCF TM模型僅需一個(gè)增量步即可得到收斂解。

圖7 L形梁結(jié)構(gòu)Fig. 7 L-shaped beam structure

圖8給出了不同拉壓模量比下，L形梁的變形后構(gòu)形圖，表1給出了采用不同技術(shù)途徑在不同拉壓模量比條件下得到的節(jié)點(diǎn)P的水平位移以及計(jì)算收斂耗時(shí)。以上仿真結(jié)果表明，即使在求解拉壓模量相差20倍的大變形問題時(shí)，本文所提出的計(jì)算框架仍能高效收斂，且能夠保證較高的精度。

圖8 變形后構(gòu)形圖Fig. 8 Configuration diagram after deformation

表1 節(jié)點(diǎn)P的水平位移與計(jì)算耗時(shí)Table 1 Horizontal displacement of nodePand calculation convergence time

book=223,ebook=127此外，如表2所示，在采用ANCF-LM模型求解不同模量大變形問題時(shí)，無法得到收斂的結(jié)果。值得注意的是，即使拉壓模量相同的條件下，ANCF-LM模型仍然無法得到收斂的結(jié)果，這是由于在相同模量時(shí)，該模型無法退化為經(jīng)典的線性彈性本構(gòu)模型。而本文提出的ANCF-TM模型求解過程中牛頓迭代收斂性良好，迭代步數(shù)和迭代耗時(shí)T隨拉壓模量比的增大略有增加。

表2 不同模型下節(jié)點(diǎn)P的求解效率Table 2 Solution efficiency of nodeP under different models

4.3 安全氣囊的充氣仿真過程

圖9所示為一方形安全氣囊的初始構(gòu)形，薄膜材料的楊氏模量為588.0 Mpa，泊松比為0.4。方膜對(duì)角線長度為L= 1 200 mm，厚度t= 1mm。考慮將氣囊的充氣過程近似為一個(gè)準(zhǔn)靜態(tài)過程，氣囊在內(nèi)部遞增的均勻壓力作用下發(fā)生變形，直至內(nèi)部壓強(qiáng)P最終增至5 000 Pa，充氣過程結(jié)束。仿真過程中，不考慮薄膜的壓縮模量。

圖9 安全氣囊的初始構(gòu)形示意圖Fig. 9 Schematic diagram of initial configuration of airbag

表3給出了基于本文所提出的計(jì)算方法在不同網(wǎng)格劃分密度下得到的安全氣囊上節(jié)點(diǎn)O、C、B的位移以及計(jì)算耗時(shí)T。結(jié)果表明，基于本文計(jì)算框架得到的位移收斂解與文獻(xiàn)[20]的研究結(jié)果基本吻合。此外，采用10×10的網(wǎng)格即可得到位移收斂解，相應(yīng)的計(jì)算耗時(shí)隨著網(wǎng)格的加密有大幅增加。

表3 不同網(wǎng)格密度下的薄膜節(jié)點(diǎn)的位移Table 3 Displacement of membrane nodes under different mesh densitied

圖10出了安全氣囊充氣完成后的應(yīng)力分布示意圖，根據(jù)應(yīng)力準(zhǔn)則[13]，對(duì)單元內(nèi)每個(gè)積分點(diǎn)的應(yīng)力狀態(tài)進(jìn)行標(biāo)記，用黑點(diǎn)標(biāo)記的薄膜區(qū)域處于褶皺狀態(tài)，用圓標(biāo)記的薄膜區(qū)域處于松弛狀態(tài)，其余未標(biāo)記的區(qū)域則處于張緊狀態(tài)，這一結(jié)果與文獻(xiàn)[21]的結(jié)論基本一致。綜上結(jié)論再次驗(yàn)證了本文所提出的計(jì)算框架在求解薄膜面外變形問題時(shí)的可行性。

圖10 安全氣囊變形后應(yīng)力分布示意圖Fig. 10 Schematic diagram of stress distribution of airbag after inflation

5 薄膜太陽帆的自旋展開動(dòng)力學(xué)模擬

作為一種新型空間可展開結(jié)構(gòu)，能夠通過自旋展開并最終穩(wěn)定的大型帆膜結(jié)構(gòu)在未來的空間探索技術(shù)應(yīng)用中具有巨大潛力。能夠準(zhǔn)確地對(duì)薄膜展開過程的動(dòng)態(tài)響應(yīng)進(jìn)行預(yù)測(cè)和模擬是薄膜太陽帆展開系統(tǒng)設(shè)計(jì)中的關(guān)鍵技術(shù)之一。這里將采用本文所提出的計(jì)算框架對(duì)薄膜太陽帆的自旋展開過程進(jìn)行分析，通過確定影響展開的關(guān)鍵因素來指導(dǎo)太陽帆的結(jié)構(gòu)設(shè)計(jì)。

初始時(shí)刻，薄膜太陽帆按照?qǐng)D11所示的折痕折疊后纏繞在與其同心的正六邊形彀輪之上，彀輪驅(qū)動(dòng)后，太陽帆隨之旋轉(zhuǎn)，并在均布于其六個(gè)角上的集中質(zhì)量塊的離心力作用下自旋展開，完全展開后的太陽帆是一個(gè)正六邊形。此時(shí)太陽帆的邊長為L= 0.635 m，厚度為t= 0.1 mm。薄膜材料的密度為ρ= 1 420 kg/m3，book=224,ebook=128楊氏模量為E= 2.5×109Pa，泊松比為v= 0.3。彀輪的邊長為d= 0.1 m，質(zhì)量為m0=1.0 kg，相對(duì)于過質(zhì)心的Z軸的慣性矩為Iz= 5.0×10–3kg m2。薄膜邊界角點(diǎn)處的集中質(zhì)量塊的質(zhì)量為m= 0.05 kg，初始時(shí)刻系統(tǒng)的角速度為w= 0.628 rad/s。

圖11 薄膜太陽帆的有限元模型Fig. 11 Finite element model of membrane solar sail

考慮到結(jié)構(gòu)的對(duì)稱性，本文僅建立1/6模型，并施加旋轉(zhuǎn)邊界條件，該有限元模型中包括35個(gè)薄膜片，薄膜片之間通過旋轉(zhuǎn)鉸連接，具體建模方法可以參考文獻(xiàn)[8]。本文采用考慮拉壓不同模量的ANCF薄膜單元和第一類拉格朗日方程對(duì)太陽帆進(jìn)行動(dòng)力學(xué)精確建模，然后利用可控?cái)?shù)值耗散的廣義-α算法對(duì)太陽帆的動(dòng)力學(xué)展開進(jìn)行求解分析。

5.1 壓縮剛度折減系數(shù)對(duì)展開過程的影響

在對(duì)太陽帆進(jìn)行仿真時(shí)考察了3種不同的剛度折減的薄膜模型（Ec=s×Et），其中s為剛度折減系數(shù)。模型一為不考慮壓縮剛度的不可壓ANCF薄膜單元（s=0）；模型二為考慮微弱抗壓剛度的ANCF薄膜單元（s=10–5）；模型三同樣為引入微弱抗壓剛度的ANCF薄膜單元（s= 10–4）。

圖12給出了采用模型一進(jìn)行仿真時(shí)薄膜太陽帆的展開過程圖，結(jié)果表明，此時(shí)的薄膜太陽帆無法完全展開并維持展開狀態(tài)，太陽帆在12 s左右展開至最大后又重新收縮至中心彀輪附近。

圖12 太陽帆的展開過程圖（s=0）Fig. 12 Deployment process diagram of solar sail （s=0）

圖13和圖14分別給出了采用模型一進(jìn)行仿真時(shí)太陽帆展開過程中的能量與角動(dòng)量演化過程。圖13表明，即使展開過程中薄膜的應(yīng)變能出現(xiàn)高頻振蕩，但是系統(tǒng)的總能量仍然是守恒的。在12 s左右時(shí)太陽帆系統(tǒng)的應(yīng)變能達(dá)到了峰值，此時(shí)太陽帆處于最大張緊狀態(tài)。圖14表明，相比于中心彀輪和集中質(zhì)量塊，薄膜在整個(gè)仿真過程中的角動(dòng)量僅發(fā)生小范圍的波動(dòng)，在仿真時(shí)間達(dá)到12 s之后，太陽帆系統(tǒng)各部件的角動(dòng)量均趨于穩(wěn)定。

圖13 太陽帆展開過程中系統(tǒng)的應(yīng)變能及總能量（s=0）Fig. 13 Strain energy and total energy of the system during solar sail deployment （s= 0）

圖14 太陽帆展開過程中系統(tǒng)各部件的角動(dòng)量（s=0）Fig. 14 Angular momentum of each component of the system during solar sail deployment （s=0）

圖15給出了采用模型二進(jìn)行仿真時(shí)薄膜太陽帆的展開過程圖，此時(shí)的薄膜太陽帆同樣無法維持最大展開狀態(tài)，太陽帆在11 s左右展開到最大半徑，之后在邊緣集中質(zhì)量區(qū)域出現(xiàn)較大幅度的面外變形。根據(jù)圖16book=225,ebook=129和圖17給出的太陽帆系統(tǒng)展開過程中的能量與和角動(dòng)量的演化過程，可以得知，系統(tǒng)的應(yīng)變能在12 s之后穩(wěn)定在一較小的范圍內(nèi)。此外，在仿真時(shí)間達(dá)到11 s左右時(shí)，太陽帆系統(tǒng)各部件的角動(dòng)量出現(xiàn)突變，之后逐漸穩(wěn)定，展開過程中系統(tǒng)的總能量和總角動(dòng)量均是守恒的。

圖15 太陽帆的展開過程圖（s=10-5）Fig. 15 Deployment process diagram of solar sail（s=10-5）

圖16 太陽帆展開過程中系統(tǒng)的應(yīng)變能及總能量（s=10-5）Fig. 16 Strain energy and total energy of the system during solar sail deployment（s=10-5）

圖17 太陽帆展開過程中系統(tǒng)各部件的角動(dòng)量（s=10-5）Fig. 17 Angular momentum of each component of the system during solar sail deployment （s=10-5）

圖18給出了采用模型三進(jìn)行仿真時(shí)薄膜太陽帆的展開過程圖，這時(shí)的薄膜太陽帆可以完全展開并維持展開狀態(tài)。太陽帆在13 s左右展開至最大半徑，之后跟隨彀輪一起穩(wěn)定旋轉(zhuǎn)，同時(shí)在邊緣集中質(zhì)量區(qū)域發(fā)生微弱的面外變形。同時(shí)根據(jù)圖19和圖20分別給出其展開過程中的能量和角動(dòng)量的演化過程可以得知，此時(shí)薄膜太陽帆系統(tǒng)的總能量和總角動(dòng)量均是守恒的。

圖18 太陽帆的展開過程圖(s=10-4)Fig. 18 Deployment process diagram of solar sail(s=10-4)

圖19 太陽帆展開過程中系統(tǒng)的應(yīng)變能及總能量（s=10-4）Fig. 19 Strain energy and total energy of the system during solar sail deployment（s=10-4）

圖20 太陽帆展開過程中系統(tǒng)各部件的角動(dòng)量（s=10-4）Fig. 20 Angular momentum of each component of the system during solar sail deployment（s=10-4）

圖21和圖22給出了基于不同壓縮剛度折減系數(shù)的薄膜模型仿真得到的太陽帆系統(tǒng)各部件的動(dòng)力學(xué)響應(yīng)。圖21表明，當(dāng)剛度折減系數(shù)為0和10–5時(shí)，太陽帆系統(tǒng)的彀輪和邊緣質(zhì)量塊無法進(jìn)入同步旋轉(zhuǎn)。而在剛度折減系數(shù)為10–4的情況下，在仿真時(shí)間達(dá)到13 s左右時(shí)彀輪和邊緣質(zhì)量塊的旋轉(zhuǎn)角度開始相等，并進(jìn)入同步旋轉(zhuǎn)。圖22表明，在考慮壓縮剛度的情況下，當(dāng)剛book=226,ebook=130度折減系數(shù)較大時(shí)，薄膜太陽帆展開至一定半徑大小之后趨于穩(wěn)定。而當(dāng)剛度折減系數(shù)較小或者不考慮壓縮剛度時(shí)，薄膜太陽帆無法維持最大展開狀態(tài)。綜上分析可知，在對(duì)薄膜太陽帆系統(tǒng)進(jìn)行展開動(dòng)力學(xué)分析時(shí)，考慮一定的壓縮剛度是必要的，這與文獻(xiàn)[4]得出的結(jié)果基本是一致的。

圖21 不同剛度折減系數(shù)時(shí)邊界質(zhì)量塊和中心彀輪的轉(zhuǎn)角Fig. 21 Rotation angle of edge mass and central hub under different stiffness reduction factors

圖22 不同剛度折減系數(shù)時(shí)的太陽帆展開半徑Fig. 22 Deployment radius of solar sail under different stiffness reduction factors

book=227,ebook=131

5.2 集中質(zhì)量和初始轉(zhuǎn)速對(duì)展開過程影響分析

為進(jìn)一步確定薄膜太陽帆的一些結(jié)構(gòu)參數(shù)對(duì)其展開動(dòng)力學(xué)特性的影響，本文又分別仿真了當(dāng)邊界集中質(zhì)量和系統(tǒng)初始轉(zhuǎn)速變化時(shí)，太陽帆系統(tǒng)各部件的動(dòng)力學(xué)響應(yīng)。

圖23和圖24分別給出了太陽帆系統(tǒng)邊界集中質(zhì)量為0.03、0.05、0.07 kg時(shí)其各部件的轉(zhuǎn)角演化過程和太陽帆展開半徑演化過程。分析得知，太陽帆分別在16、13、11 s左右展開至最大半徑，同時(shí)中心彀輪和邊界質(zhì)量塊的旋轉(zhuǎn)角度開始趨于同步，展開過程的構(gòu)形演化圖類似圖18。此外，綜合分析得知，邊界集中質(zhì)量對(duì)太陽帆系統(tǒng)的最終展開半徑大小以及完全展開后發(fā)生的面外變形均會(huì)產(chǎn)生一定影響。

圖23 不同集中質(zhì)量時(shí)邊界質(zhì)量塊和中心彀輪的轉(zhuǎn)角Fig. 23 Rotation angle of edge mass and central hub under different concentrated masses

圖24 不同集中質(zhì)量時(shí)的太陽帆展開半徑Fig. 24 Deployment radius of solar sail with different concentrated masses

圖25和圖26給出了系統(tǒng)初始轉(zhuǎn)速分別為0.471、0.628、0.785 rad/s時(shí)太陽帆展開過程中各部件的動(dòng)力學(xué)響應(yīng)。結(jié)果表明，太陽帆分別在16.5、13、10 s左右展開至最大半徑，同時(shí)中心彀輪和邊界質(zhì)量塊開始進(jìn)入同步旋轉(zhuǎn)，其展開過程的構(gòu)形演化類似圖18。綜合分析得知，系統(tǒng)初始轉(zhuǎn)速對(duì)太陽帆展開至最大半徑所耗費(fèi)的時(shí)間以及完全展開后發(fā)生的面外變形均會(huì)產(chǎn)生顯著影響。

圖25 不同初始轉(zhuǎn)速時(shí)邊界質(zhì)量塊和中心彀輪的轉(zhuǎn)角Fig. 25 Rotation angle of edge mass and central hub at different initial speeds

圖26 不同初始轉(zhuǎn)速時(shí)的太陽帆展開半徑Fig. 26 Deployment radius of solar sail at different initial speeds

6 結(jié) 論

本文整合ANCF建模方法與張力場(chǎng)理論，針對(duì)薄膜結(jié)構(gòu)拉壓不同模量的宏觀特性，提出了一種高效的建模與計(jì)算方法，為空間帆膜結(jié)構(gòu)的展開動(dòng)力學(xué)模擬提供了有效途徑。

1）針對(duì)薄膜太陽帆系統(tǒng)展開過程中顯著的非線性動(dòng)力學(xué)特性，基于張力場(chǎng)理論與不同模量彈性理論，推導(dǎo)了單元精確的切向剛度矩陣，提出了一種考慮不同模量的ANCF薄膜單元計(jì)算框架。其中考慮了薄膜單元應(yīng)力場(chǎng)對(duì)其本構(gòu)關(guān)系的影響，準(zhǔn)確地描述了薄膜結(jié)構(gòu)在拉壓過渡區(qū)的非線性特性，有效提高了單元在牛頓迭代過程中的收斂性。

2）基于上述框架，采用第一類拉格朗日方程建立了六邊形薄膜太陽帆的自旋展開動(dòng)力學(xué)模型，并通過可控?cái)?shù)值耗散的廣義-α算法求解運(yùn)動(dòng)方程。根據(jù)數(shù)值模擬結(jié)果，分析了太陽帆不同結(jié)構(gòu)設(shè)計(jì)參數(shù)對(duì)其展開book=228,ebook=132動(dòng)力學(xué)特性的影響，為太陽帆類空間帆膜結(jié)構(gòu)的系統(tǒng)設(shè)計(jì)和折展動(dòng)力學(xué)仿真提供了一定的理論指導(dǎo)。

3）未來值得進(jìn)一步研究的問題和方向包括：薄膜太陽帆完全展開之后，帆面產(chǎn)生的微弱面外變形以及振蕩過程仍未得到有效解決；現(xiàn)有的針對(duì)薄膜太陽帆自旋展開的動(dòng)力學(xué)建模與仿真方法在計(jì)算效率和分析精度方面仍然存在不足，隨著薄膜單元網(wǎng)格的加密，難以針對(duì)大型帆膜結(jié)構(gòu)進(jìn)行高效模擬。