航天器中頻段力學環(huán)境預(yù)示的混合線連接建模方法

朱衛(wèi)紅，馬興瑞，韓增堯，鄒元杰

(中國空間技術(shù)研究院，北京100094)

0 引言

準確預(yù)示航天器的力學環(huán)境是指導(dǎo)航天器設(shè)計和地面試驗方案與條件制定的重要依據(jù)，但是航天器在發(fā)射過程中經(jīng)歷的力學環(huán)境頻率寬，一般在10 Hz～10 kHz，很難采用單一的分析方法實現(xiàn)精確的力學環(huán)境預(yù)示［1］。在低頻段航天器結(jié)構(gòu)和聲腔的模態(tài)稀疏，有限元和邊界元等基于單元離散技術(shù)的方法最為常用;而在高頻段，結(jié)構(gòu)和聲腔的模態(tài)密集且隨機特性影響突出，因此，統(tǒng)計能量分析和能量有限元分析等都是比較有效的預(yù)示手段。但是在低頻和高頻中間存在一個過渡頻段(一般對于航天器系統(tǒng)級模型大約在300 Hz～1 kHz)，在該頻段航天器各子結(jié)構(gòu)或部組件之間的動力學特性差異大，此時傳統(tǒng)低頻方法和高頻方法都無法準確地進行預(yù)示，這就是航天器中頻段的力學環(huán)境預(yù)示問題。目前中頻段的響應(yīng)預(yù)示是航天器結(jié)構(gòu)力學環(huán)境預(yù)示的難點，也是國內(nèi)外研究的熱點［2］，因此展開航天器中頻段的預(yù)示方法研究具有重要的理論價值和工程應(yīng)用價值。

為了解決中頻段的力學環(huán)境預(yù)示問題，2005年Shorter等［3］提出了一種混合FE-SEA中頻預(yù)示方法(Hybrid Finite Element-Statistic Energy Method)。該方法首先將系統(tǒng)根據(jù)動力學特性劃分為確定性子系統(tǒng)和隨機子系統(tǒng)，然后確定性子系統(tǒng)采用有限元方法建模，隨機子系統(tǒng)采用SEA方法建模，最后耦合求解得到系統(tǒng)不同子系統(tǒng)的響應(yīng)。該方法提出后受到了廣泛關(guān)注，目前國外已將其廣泛應(yīng)用于航空航天領(lǐng)域［4－6］。我國航天工業(yè)部門也在“十一五”初期啟動了該方法的理論和應(yīng)用研究，目前已完成了該方法的基本理論推導(dǎo)和部分仿真分析及試驗驗證［7－9］，但是結(jié)合目前我國航天工程中頻預(yù)示的實際需求，仍然需要針對該方法從建模理論、工程應(yīng)用和試驗驗證等方面展開更深入的研究［2］。

對于大型復(fù)雜航天器，不同部組件的連接形式通常復(fù)雜多樣，總體可歸納為點、線、面三種基本連接形式，其中線連接是航天器結(jié)構(gòu)中一種比較普遍的連接形式，如星體壁板的交界處，又如衛(wèi)星與火箭連接的星箭界面等。如果線連接的兩類結(jié)構(gòu)剛度差異大，這就是典型的中頻預(yù)示問題。根據(jù)混合FESEA理論，需采用有限元方法對剛性較大的部件建模，而柔性組件則采用統(tǒng)計能量方法建模，然后通過混合線連接耦合進行求解分析。目前國內(nèi)外公開可查閱的混合線連接建模的資料非常少，因此本文針對航天器結(jié)構(gòu)中頻預(yù)示的混合線連接建模理論展開研究，最后以典型結(jié)構(gòu)為例驗證了本文方法的有效性。

1 基于混合線連接模型的FE-SEA分析

混合FE-SEA基本理論認為統(tǒng)計能量子系統(tǒng)的能量由直接場和混響場兩部分能量組成。直接場能量為通過連接輸入到子系統(tǒng)后未經(jīng)過任何邊界反射的能量，混響場能量為經(jīng)過多次反射后的能量，圖1為線連接結(jié)構(gòu)的直接場和混響場示意圖。

假設(shè)統(tǒng)計能量子系統(tǒng)的混合線連接m的位移自由度為q(m)，則混合線連接的動力學響應(yīng)可表示為:

圖1 線連接的直接場和混響場Fig.1 Direct field and diffuse reverberant field for line junction

上標－H為共軛轉(zhuǎn)置求逆運算符，力互譜矩陣Sff包括外載荷互譜和混響受擋力互譜:

式中:Im(·)為取虛部，E為統(tǒng)計能量子系統(tǒng)的混響能量，ω為分析頻率，n為模態(tài)密度。方程描述了統(tǒng)計能量子系統(tǒng)的能量響應(yīng)與其在混合線連接m處產(chǎn)生的反力之間的關(guān)系。同時基于該關(guān)系可建立耦合后的統(tǒng)計能量模型的功率流平衡方程［11］:

式中:式中:Ddir，j為第j個統(tǒng)計能量子系統(tǒng)所有混合線連接直接場動剛度矩陣之和，ηd，j為有限元模型對統(tǒng)計能量模型j的等效損耗因子，ηje為SEA子系統(tǒng)j和SEA子系統(tǒng)e間的等效耦合損耗因子，Pjin為施加在統(tǒng)計載荷輸入到SEA子系統(tǒng)j的等效輸入功率，nj為SEA子系統(tǒng)j的模態(tài)密度。分析時首先求解耦合后的功率流平衡方程(6)，求解得到統(tǒng)計能量模型的能量響應(yīng)，然后由式(5)得混合連接邊界處的混響受擋力，代入方程(3)求得有限元模型的響應(yīng)。

從上述求解流程可以看出，混合FE-SEA分析關(guān)鍵在于混合連接的建模，即混合線連接的直接場動剛度矩陣。下面將針對如何建立混合線連接的直接場動剛度矩陣展開研究。

2 混合線連接建模

2. 1 波數(shù)空間下線連接動剛度矩陣

航天器結(jié)構(gòu)以具有大面質(zhì)比的板殼結(jié)構(gòu)為主，如星體壁板結(jié)構(gòu)和太陽翼等，因此本文重點研究板殼結(jié)構(gòu)的混合線連接模型。假設(shè)板j通過線連接與其他結(jié)構(gòu)相連，定義其局部坐標系為oxjyjzj，如圖2所示。線連接沿xj軸，其運動可由沿三個方向的平動［uj，vj，wj］T與繞 xj軸的轉(zhuǎn)角 θj描述，這里主要考慮面外運動［wj，θj］T，面內(nèi)運動的建模將不再贅述。根據(jù)各向同性薄板理論，板的橫向自由運動的控制方程為:

圖2 線連接的局部坐標系Fig.2 The local coordinate of line junction

波從線連接處向板結(jié)構(gòu)內(nèi)部傳遞能量，波的形式可寫為exp(－i kxj+μyj+iωt)，負號表示波沿xj軸正方向傳播，k為該方向波數(shù)，則橫向位移可寫為:

式中:α為橫向位移的幅值，將式(8)代入式(7)則可得出:

式中:板的彎曲波數(shù) kB=(ρjhjω2/Dj)1/4。由方程可求出4個μ值，由于yj方向的波只能沿正方向傳播，因此μ只能取負值。當k≥kB時，μ取兩個負實數(shù)，這表明此時yj方向的波為近場波，不傳遞能量;當k＜kB時，μ取一個負實數(shù)和一個負虛數(shù)，此時yj方向的波為一個近場波和一個遠場波。假設(shè)μ的兩個可能值分別為μB1和μB2，相對應(yīng)的兩類波的橫向位移幅值分別為α1和α2，則板的橫向位移可以表示為:

xj方向的轉(zhuǎn)角θj可由橫向位移的偏導(dǎo)給出:

考慮到線連接位于xj軸，省略時間項，則由式(10)和式(11)可求出:

線連接上xj軸的力矩與橫向力可表示為:

式中:v0為板的泊松比。將式(10)～(12)代入式(13)(14)中即可求得波數(shù)空間下線連接的動剛度矩陣為［12］:

式(15)即為應(yīng)用波動理論得到的波數(shù)空間下線連接處的動剛度矩陣。注意式(15)為基于波數(shù)坐標k建立的，因此無法應(yīng)用到混合FE-SEA分析中。另一方面從推導(dǎo)過程中可以看出，在基于波動理論推導(dǎo)線連接動剛度矩陣時，線連接沿xj方向為無限長，所以需要進一步處理。

2. 2 混合線連接直接場動剛度矩陣

由于有限元模型是基于節(jié)點坐標系的，所以需要將波數(shù)空間下的動剛度矩陣變換到節(jié)點坐標系下。設(shè)基于廣義坐標an的線連接的面外運動可表示為:

式中:U(x)=［w(x)，θ(x)］T，φn(x)線連接的面外位移與轉(zhuǎn)角的形函數(shù)向量。假設(shè)在線連接外的位移為0，即線連接具有固支邊界條件，則線連接上直接場動剛度矩陣的第mn項{Ddir}mn可表示為［13］:

式中:L為實際結(jié)構(gòu)中線連接的長度，fn(x)為線連接處的外力 f(x)在 φn(x)上的分量，對形函數(shù)φm(x)應(yīng)用空間傅立葉變換:

在波數(shù)空間下線連接處位移形函數(shù)與其引起的相應(yīng)分布力之間的關(guān)系可表示為:

式中:D∞(k)為波數(shù)空間下半無限板線連接的動剛度陣，由式(15)給出。對式(19)應(yīng)用傅立葉逆變換可得:

將式(20)代入式(17)有:

式(21)給出了線連接直接場動剛度矩陣的第mn項的計算格式，可以看出線連接的直接場動剛度矩陣與形函數(shù)相關(guān)，形函數(shù)的構(gòu)造方法是本文研究的重點。

3 形函數(shù)構(gòu)造及變換

由式(21)可以看出，形函數(shù)是直接場動剛度矩陣的計算的關(guān)鍵參數(shù)。本文提出了三種位移形函數(shù)構(gòu)造方法，下面將對三種方法逐一進行介紹。

3. 1 模態(tài)線性插值形函數(shù)

形函數(shù)的選取可基于線連接處的有限元模態(tài)分析結(jié)果。首先對混合模型中的有限元模型進行模態(tài)分析，提取線連接上與其節(jié)點相關(guān)的N階模態(tài)，然后對模態(tài)進行線性插值獲得第m個形函數(shù)φm(x):

式中:系數(shù)am和bm通過對第m個模態(tài)進行線性插值獲得，對式(22)進行傅里葉變換:

式(23)即為由模態(tài)j線性插值構(gòu)造得到波數(shù)空間下的形函數(shù)。可以看出模態(tài)插值形函數(shù)的計算比較復(fù)雜，首先必須進行有限元分析提取線連接節(jié)點的模態(tài)，然后進行線性插值后求其傅里葉變換，最后才能代入式(21)中求解節(jié)點坐標系的動力學剛度矩陣。

3. 2 三角波線性插值形函數(shù)

模態(tài)函數(shù)的插值和傅里葉變換計算繁瑣復(fù)雜，形函數(shù)的數(shù)目取決于模態(tài)的數(shù)目，因此當模型模態(tài)密集且分析頻率范圍非常寬時，采用該方法建立線連接模型就非常困難了。為了解決這個問題，本節(jié)構(gòu)造一種三角波線性插值形函數(shù)。構(gòu)造方法如下:第m個形函數(shù)定義為線連接上節(jié)點m處的位移為1，其余節(jié)點位移為0，節(jié)點間的位移采用線性插值，如圖3所示。

圖3 三角波形函數(shù)Fig.3 Triangle wave shape function

若節(jié)點的位置為x=0時，形函數(shù)可寫為:

式中:Δl為節(jié)點間的距離。該三角波形函數(shù)的空間傅里葉變換可寫為:

根據(jù)傅里葉變換的時移特性，若第m個形函數(shù)的位置為x=xm，則有:

式(26)給出了三角波形函數(shù)的表達式。從推導(dǎo)過程可以看出，基于三角波插值的形函數(shù)特點在于形函數(shù)的數(shù)目由線連接上的節(jié)點數(shù)目決定，這是其優(yōu)于模態(tài)插值形函數(shù)的一個重要特點，尤其是在寬頻帶分析時。但是注意到式(21)是奇異積分，因此采用上述兩種形函數(shù)構(gòu)造方法可能存在積分不收斂的問題。

3. 3 基于Shannon小波的形函數(shù)

為了解決式積分由可能存在不收斂的問題，本節(jié)提出了一種基于Shannon小波 (Sinc函數(shù))的線連接位移形函數(shù)構(gòu)造方法，該函數(shù)的傅里葉變換在波數(shù)空間下是窗函數(shù)(如圖4所示)，可有效簡化線連接建模并解決積分奇異問題。

圖4 Shannon小波函數(shù)及其傅里葉變換Fig.4 Shannon wavelet and it＇s Fourier transform

定義以Shannon小波函數(shù)為基的形函數(shù)為:

式中:ks與節(jié)點相關(guān)，本文取ks=2π/Δl。對上式進行傅里葉變換有:

由傅里葉變換的時移定理，第m個形函數(shù)φm(x)的傅里葉變換為:

將式(28)和式(29)代入式(21)可以得到基于Shannon小波形函數(shù)的線連接動剛度矩陣為:

式(30)表明利用Shannon小波函數(shù)的傅里葉變換特性，直接場的動剛度矩陣的積分變?yōu)橛薪绶e分，且積分只與節(jié)點間的相對位置有關(guān)，因此簡化了計算，可采用常規(guī)的積分算法求解式(30)。

4 仿真試驗

為了校驗本文提出的混合線連接模型的有效性，設(shè)計如圖5所示的星體壁板與連接結(jié)構(gòu)的組合模型對本文提出的混合線連接模型進行校驗。模型為星體壁板連接的簡化模型，兩塊星體壁板結(jié)構(gòu)通過復(fù)雜的剛性結(jié)構(gòu)連接，每個壁板結(jié)構(gòu)與剛性連接結(jié)構(gòu)之間通過線連接進行耦合。星體壁板結(jié)構(gòu)楊氏模量為71 GPa，密度2700 kg/m3，泊松比 0.3296，剛性連接結(jié)構(gòu)楊氏模量為 210 GPa，密度為 7800 kg/m3，泊松比0.3125，整個結(jié)構(gòu)的損耗因子為1%，模型的幾何參數(shù)如表1所示。載荷為點載荷，垂直施加在其中一個星體壁板結(jié)構(gòu)上，分析頻率范圍為1 Hz～1 kHz?？紤]到整個系統(tǒng)的動力學特性，柔性較大的板面外位移采用SEA建模，剛性較大的板的面內(nèi)位移和連接結(jié)構(gòu)采用FEM建模，兩種模型通過兩個混合線連接耦合，如圖6所示。參考算法采用基于有限元的能量流分析與蒙特卡洛仿真(Monte Carlo Simulation)相結(jié)合的方法［14］，并與國外商用軟件 VA One［15］的計算結(jié)果進行對比。

圖5 組合結(jié)構(gòu)模型Fig.5 Combined structure model

圖6 混合線連接模型Fig.6 Hybrid line junction model

表1 結(jié)構(gòu)的幾何參數(shù)Table 1 Geometrical parameters of the structure

為了描述中高頻響應(yīng)對參數(shù)攝動的敏感性，在壁板結(jié)構(gòu)中引入不確定性參數(shù)，本文中在每塊板上附加20個隨機分布的集中質(zhì)量塊，集中質(zhì)量塊總質(zhì)量為壁板結(jié)構(gòu)的15%，以有限元能量流分析—蒙特卡洛仿真得到系統(tǒng)的平均響應(yīng)作為參考。圖7為兩個統(tǒng)計能量子系統(tǒng)間的耦合損耗因子計算結(jié)果。仿真結(jié)果表明，本文提出的基于Shannon小波的線連接建模方法和VA One的計算結(jié)果均與Monte Carlo仿真的結(jié)果都相吻合。但是VA One計算結(jié)果偏于保守，尤其是在非共振峰位置。

圖8給出了壁板子系統(tǒng)的無量綱響應(yīng)結(jié)果(子系統(tǒng)能量與系統(tǒng)輸入能量的比值)。從圖中可以看出，對于能量輸入的壁板結(jié)構(gòu)1，本文方法和 VA One均與Monte Carlo試驗相吻合，但是對壁板結(jié)構(gòu)2，本文提出的方法在低頻處與Monte Carlo試驗存在誤差，其余頻段都比較吻合，而VA One誤差在非共振峰位置相對較大。低頻處的差異主要是由于在混合線連接建模中引入的固支邊界假設(shè)造成的。

圖7 基于Shannon小波方法的耦合損耗因子Fig.7 The coupling loss factor based on the Shannon wavelet method

圖8 基于Shannon小波方法的壁板響應(yīng)Fig.8 The plates response based on the Shannon wavelet method

圖9 不同形函數(shù)的壁板之間的耦合損耗因子計算結(jié)果Fig.9 The coupling loss factor of two plates with different shape function

圖10 不同形函數(shù)的壁板響應(yīng)預(yù)示結(jié)果對比Fig.10 The plates response results with different shape function

圖9 和圖10為基于不同形函數(shù)的線連接模型計算得到的兩個壁板之間的耦合損耗因子和壁板結(jié)構(gòu)響應(yīng)的對比曲線。從圖中可以看出，三種形函數(shù)模型下得到的結(jié)果基本吻合，但是模態(tài)形函數(shù)的誤差相對較大，這是由于模態(tài)截斷和線性插值引入的誤差導(dǎo)致的。從結(jié)果的對比可以看出，本文提出的三種形函數(shù)都可以有效建立混合線連接模型，對結(jié)構(gòu)的響應(yīng)進行正確預(yù)示。但是基于Shannon小波形函數(shù)的線連接模型構(gòu)造最為簡單，能夠有效簡化計算過程，解決奇異積分不收斂的問題，因此更適合于航天器工程應(yīng)用。

5 結(jié)論

本文研究了航天器結(jié)構(gòu)中頻段力學環(huán)境預(yù)示的混合線連接建模方法，并針對典型結(jié)構(gòu)展開了仿真校驗，主要內(nèi)容和結(jié)論如下:

(1)采用傅里葉變換方法，將波數(shù)空間下的線連接動剛度矩陣變換為有限元節(jié)點坐標系下的線連接直接場動剛度矩陣，傅里葉變換技術(shù)同樣也實現(xiàn)了將無限線連接模型轉(zhuǎn)化為符合實際結(jié)構(gòu)的有限線連接模型。

(2)提出的三種形函數(shù)構(gòu)造方法中，Shannon小波形函數(shù)線連接模型具有建模簡單、可解決動剛度矩陣積分奇異的問題等特點;模態(tài)插值線連接模型建模相對比較復(fù)雜，特別是當分析頻帶非常寬且模態(tài)密集時，該模型的求解效率急劇下降。

(3)算例仿真分析表明，本文提出的混合線連接建模方法具有較好的預(yù)示精度。同時對比三種形函數(shù)線連接模型的預(yù)示結(jié)果可以看出預(yù)示結(jié)果基本吻合，其中模態(tài)線性插值函數(shù)誤差相對較大。

(4)本文的驗證對象是簡化的星體壁板連接結(jié)構(gòu)，下一步工作應(yīng)當從算法優(yōu)化、更加復(fù)雜的航天器結(jié)構(gòu)的應(yīng)用和試驗驗證兩方面展開深入研究。

［1］ 馬興瑞，于登云，韓增堯，等.星箭力學環(huán)境分析與試驗技術(shù)研究進展［J］.宇航學報，2006，27(3):323－331.［Ma Xing-rui， Yu Deng-yun， Han Zeng-yao， et al. Research evolution on the satellite-rocket mechanical environment analysis＆ test technology［J］.Journal of Astronautics，2006，27(3):323 －331.］

［2］ 馬興瑞，韓增堯，鄒元杰，等.航天器力學環(huán)境分析與條件設(shè)計研究進展［J］.宇航學報，2012，33(1):1－12.［Ma Xing-rui，Han Zeng-yao，Zou Yuan-jie，et al. Review and assessment of spacecraft mechanical environment analysis and specification determination［J］.Journal of Astronautics，2012，33(1):1 －12.］

［3］ Shorter P J，Langley R S.Vibro-acoustic analysis of complex systems［J］.Journal of Sound and Vibration，2005，288(3):669－699.

［4］ Jeffrey M L，Cotoni V.Vibro-acoustic response of the NASA ACTS spacecraft antenna to launch acoustic excitation［R］.Washington D C，Glenn Research Center，NASA/TM －2008－215168，1－15，2008.

［5］ Chen SM，Wang D F，Zan JM.Interior noise prediction of the automobile based on hybrid FE-SEA method［J］.Mathematical Problems in Engineering，2011，2011(4):1－20.

［6］ Prasanth S，Charpentier A，F(xiàn)ukui K.Using the hybrid FE-SEA model of a trimmed full vehicle to reduce structure born noise from 200Hz to 1kHz［C］.12th Symposium on International Automotive Technology，Pune，India，2011.

［7］ 張瑾.FE-SEA方法在航天器力學環(huán)境預(yù)示中的應(yīng)用研究［D］.北京:中國空間技術(shù)研究院，2011.［Zhang Jin.Spacecraft mid-frequency dynamic predication using FE-SEA hybrid method［D］. Beijing:China Academy of Space Technology，2011.］

［8］ 鄒元杰，韓增堯，張瑾.航天器全頻域力學環(huán)境預(yù)示技術(shù)研究進展［J］.力學進展，2012，42(4):445－454.［Zou Yuanjie，Han Zeng-yao，Zhang Jin.Research progress on fullfrequency predication techniques of spacecraft’s mechanical environment［J］.Advances in Mechanics，2012，42(4):445 －454.］

［9］ 朱衛(wèi)紅，馬興瑞，韓增堯.航天器中頻力學環(huán)境預(yù)示技術(shù)研究進展［J］.航天器工程，2014，23(1):110－117.［Zhu Wei-hong，Ma Xing-rui，Han Zeng-yao.Research evolution on mid-frequency mechanical environment predication of spacecraft［J］.Spacecraft Engineering，2014，23(1):110 －117.］

［10］ Shorter PJ，Langley R S.On the reciprocity relationship between direct field radiation and diffuse reverberant loading［J］.The Journal of the Acoustic Society of America，2005，117(1):85－95.

［11］ Langley R S，Cordioli J A.Hybrid deterministic– statistical analysis of vibro-acoustic systems with domain coupling on statistical components［J］.Journal of Sound and Vibration，2009，321(1):893－912.

［12］ Langley R S，Heron K H.Elastic wave transmission through plate/beam junction［J］.Journal of Sound and Vibration，1990，143(2):241－253.

［13］ Cotoni V，Shorter P J.Numerical and experimental validation of a hybrid finite element-statistical energy analysis method［J］.The Journal of the Acoustic Society of America，2007，122(1):259－270.

［14］ Mace B R，Shorter P J.Energy flow models from finite element analysis［J］.Journal of Sound and Vibration，2000，233(3):369－389.

［15］ VA One 2010.5 user’s guide［Z］.ESI Group，Paris，F(xiàn)rance，2010.