廖云龍,吳佩沛,尤 罡
(1. 上??臻g推進(jìn)研究所,上海 201112; 2. 上??臻g發(fā)動機(jī)工程技術(shù)研究中心,上海 201112)
空間推進(jìn)系統(tǒng)設(shè)計領(lǐng)域廣泛使用具有螺紋連接、箍帶連接的結(jié)構(gòu),易出現(xiàn)非線性模態(tài)問題,如推進(jìn)系統(tǒng)用電磁閥、氣瓶、管路組件會因螺紋、箍帶松動造成模態(tài)頻率降低。航天產(chǎn)品對力學(xué)環(huán)境條件較為敏感,而非線性模態(tài)缺乏成熟的仿真方法,這一矛盾使得非線性模態(tài)機(jī)理研究愈發(fā)迫切。如何進(jìn)行非線性模態(tài)仿真,成為空間推進(jìn)系統(tǒng)設(shè)計過程中亟待解決的問題。
非線性模態(tài)研究發(fā)展過程中,Rosenberg的工作具有重要意義[1-2]。Rosenberg 于20世紀(jì)60年代初,在對非線性系統(tǒng)的自由振動進(jìn)行研究時引入了非線性模態(tài)的概念。20世紀(jì)90年代,Shaw等提出了阻尼非線性模態(tài),將相空間的不變流形定義為非線性模態(tài),系統(tǒng)初始狀態(tài)符合某一非線性模態(tài),隨后按照該模態(tài)運(yùn)動[3-4]。劉鍊生等推廣了相似模態(tài)的概念[5-6],使其能夠應(yīng)用于非線性系統(tǒng)。李誠等提出基于譜單元的伽遼金法求解不變流形定義下的非線性模態(tài)曲面,獲得兩自由度非線性系統(tǒng)更高精度的非線性模態(tài)曲面解[7]。胡曉君應(yīng)用混合攝動伽遼金法進(jìn)行單自由度體系的非線性隨機(jī)結(jié)構(gòu)自由振動數(shù)值分析[8]。近20年來,非線性模態(tài)理論發(fā)展仍停留在簡單自由度系統(tǒng)分析階段。
非線性模態(tài)分析由于其復(fù)雜性,使其長期處于理論研究狀態(tài),工程實際多采用模態(tài)試驗與線性模態(tài)仿真結(jié)合的方法。代鵬等進(jìn)行了脈沖等離子推進(jìn)系統(tǒng)模態(tài)分析和試驗驗證[9]。杜大華等進(jìn)行了火箭發(fā)動機(jī)渦輪盤模態(tài)影響因素與振動安全性分析[10]。田彤輝等研究了沖擊荷載下級間螺栓法蘭連接結(jié)構(gòu)失效實驗與數(shù)值仿真[11]。
直接應(yīng)用非線性理論進(jìn)行模態(tài)仿真存在困難,但模擬試驗過程通過功率譜對系統(tǒng)動態(tài)參數(shù)進(jìn)行識別,將瞬態(tài)仿真方法與信號分析方法相結(jié)合,可能成為解決非線性模態(tài)問題的另一種途徑。功率譜密度(power spectral density,PSD)分析法應(yīng)用于系統(tǒng)參數(shù)識別始于文獻(xiàn)[12]。在時變系統(tǒng)或非線性系統(tǒng)參數(shù)識別相關(guān)研究領(lǐng)域,張潔等進(jìn)行了機(jī)械變速過程的瞬時參數(shù)識別及頻譜分析[13];王豪等應(yīng)用圖像處理技術(shù),提出一種自適應(yīng)時頻脊線的方法,對時變結(jié)構(gòu)在隨機(jī)激勵下的瞬時頻率進(jìn)行了辨識[14];王曉敏等研究了在白噪聲激勵下非線性聲學(xué)超材料的隨機(jī)響應(yīng)[15];許文峰等研究了具有非線性剛度邊界的桿梁結(jié)構(gòu)的動態(tài)特性[16]。傳統(tǒng)研究認(rèn)為瞬態(tài)時域仿真應(yīng)用存在計算量大、求解過程慢、數(shù)據(jù)處理過程復(fù)雜等問題,但近年來隨著計算機(jī)能力的快速提升,基于瞬態(tài)時域法與信號后處理進(jìn)行動力學(xué)仿真在多場景下已具備工程應(yīng)用條件。聶肇坤等對運(yùn)載火箭艙段連接結(jié)構(gòu)進(jìn)行了簡化建模,采用沖擊、周期信號時域非線性仿真,對時域結(jié)果進(jìn)行FFT變換后與精細(xì)模型進(jìn)行了對比驗證,使得簡化模型在保留線性計算效率的同時接近精細(xì)模型的計算精度[17]。樂晨等進(jìn)行了運(yùn)載火箭大直徑螺栓承受拉彎耦合載荷的失效判據(jù)研究[18]。
基于計算機(jī)能力的提升與功率譜在模態(tài)試驗中應(yīng)用的優(yōu)勢,本文提出一種全新的瞬態(tài)動力學(xué)仿真與功率譜密度分析結(jié)合的方法,應(yīng)用于空間推進(jìn)系統(tǒng)電磁閥線圈螺紋松動非線性模態(tài)仿真。在時域進(jìn)行接觸非線性仿真,并將時域結(jié)果變換為功率譜,在頻域進(jìn)行非線性模態(tài)分析,關(guān)注重點(diǎn)為接觸非線性對模態(tài)頻率的影響。
某型號電磁閥是推進(jìn)系統(tǒng)關(guān)鍵組件,線圈元件通過外套螺母連接于電磁閥本體,飛行試驗過程中因隨機(jī)振動導(dǎo)致線圈松動,出現(xiàn)模態(tài)偏移,本文以該閥門為研究對象,采用功率譜分析法進(jìn)行非線性模態(tài)仿真方法研究。閥門隨機(jī)試驗條件以PSD形式給出,如表1所示,隨機(jī)試驗狀態(tài)如圖1所示。
圖1 電磁閥振動試驗Fig.1 Solenoid valve vibration test
表1 隨機(jī)試驗條件Tab.1 Random test conditions
采用隨機(jī)試驗條件作為結(jié)構(gòu)的激勵條件,可在較短的計算時間內(nèi)激發(fā)結(jié)構(gòu)多階模態(tài)響應(yīng),相較正弦掃描等試驗條件更適用于時域仿真分析。但應(yīng)用隨機(jī)試驗條件進(jìn)行時域仿真需要對PSD進(jìn)行時域擬合,擬合出的隨機(jī)試驗條件具備相同的統(tǒng)計性質(zhì),但并不唯一。隨機(jī)瞬態(tài)動力學(xué)分析首先需要進(jìn)行時域隨機(jī)數(shù)據(jù)擬合方法研究。
在現(xiàn)代信號分析中,PSD分析是研究具有時不變、各態(tài)歷經(jīng)、平穩(wěn)隨機(jī)信號的重要途徑,PSD可反映隨機(jī)信號各頻率成分功率的分布情況,反映各階模態(tài)的共振效果。隨機(jī)信號的研究往往需要轉(zhuǎn)換到頻域進(jìn)行,而理論隨機(jī)信號通常是無限長的,需要使用有限長信號對真實PSD進(jìn)行估計。經(jīng)典譜估計方法包括維納-辛欽定理、周期圖法、Welch法等,其中周期圖法為應(yīng)用廣泛且最具代表性的方法,其基本原理是對觀測到的數(shù)據(jù)進(jìn)行快速傅里葉變換(fast Fourier transform,FFT),然后取模的平方作為真實PSD估計。
取平穩(wěn)隨機(jī)信號X(n)的有限個觀察點(diǎn)X(0),X(1),…,X(n),則傅里葉變換過程為
(1)
周期圖法功率譜為
(2)
瞬態(tài)動力學(xué)仿真計算量較大,仿真時間較長,應(yīng)在滿足仿真需求的前提下盡量減少仿真資源的使用,減少仿真時間。隨機(jī)條件最低頻率為10 Hz,保證10倍最低頻率采樣,擬合總時長1 s。
對隨機(jī)試驗條件進(jìn)行時域數(shù)據(jù)擬合與驗證過程如圖2所示。
圖2 時域隨機(jī)條件合成流程Fig.2 Time domain random conditional synthesis process
依據(jù)圖2(a)流程進(jìn)行擬合驗證,首先由計算機(jī)生成兩列1 s均勻分布的白噪聲UN、VN,應(yīng)用文獻(xiàn)[19]的方法獲取正態(tài)分布的數(shù)據(jù)列XN,即
(3)
將數(shù)據(jù)通過截止頻率為隨機(jī)振動條件最高頻率10倍的低通濾波器;對XN進(jìn)行FFT變換后在頻域進(jìn)行增益調(diào)整,并采用周期圖法進(jìn)行PSD計算,使其符合試驗曲線,再進(jìn)行逆FFT變換,擬合的時域信號如圖2(b)所示。
圖2(c)為擬合時域加速度PSD曲線與標(biāo)準(zhǔn)PSD曲線對比。對擬合時域數(shù)據(jù)進(jìn)行概率密度統(tǒng)計[見圖2(d)],表明擬合數(shù)據(jù)具備較好的正態(tài)分布特征。在時域?qū)M合數(shù)據(jù)進(jìn)行數(shù)學(xué)期望與均方根計算,如表2所示。其中時域數(shù)據(jù)數(shù)學(xué)期望接近0,均方根為12.2g,與表 1試驗條件均方根加速度符合較好,作為后續(xù)瞬態(tài)動力學(xué)仿真的輸入。
表2 數(shù)據(jù)統(tǒng)計性質(zhì)Tab.2 Statistical properties of data
為分析應(yīng)用瞬態(tài)動力學(xué)方法進(jìn)行非線性模態(tài)分析的有效性,同時進(jìn)行非線性瞬態(tài)動力學(xué)模態(tài)仿真與基于線性模態(tài)疊加法的隨機(jī)振動仿真。由于計算原理不同,仿真將依據(jù)不同的流程,如圖3所示。瞬態(tài)動力學(xué)仿真過程使用擬合的隨機(jī)時域條件作為輸入,并將結(jié)果轉(zhuǎn)換為PSD曲線;而線性隨機(jī)仿真首先進(jìn)行線性模態(tài)仿真,然后直接應(yīng)用模態(tài)疊加法計算PSD曲線。
圖3 PSD仿真對比流程Fig.3 PSD simulation and comparison process
電磁閥線圈出現(xiàn)非線性模態(tài)的主要原因在于固定螺紋松動,能否在仿真方法中考慮螺紋松動過程出現(xiàn)的接觸非線性,從仿真原理上考慮接觸對模態(tài)結(jié)果的影響,是決定非線性模態(tài)仿真方法是否可行的關(guān)鍵。
基于瞬態(tài)顯式動力學(xué)分析結(jié)果進(jìn)行非線性模態(tài)仿真,在時間維度上進(jìn)行迭代求解。瞬態(tài)顯式動力學(xué)計算原理為應(yīng)用中心差分法對物體運(yùn)動過程進(jìn)行積分,其迭代求解的基本方程如式(4)與式(5)所示。
(4)
(5)
(6)
瞬態(tài)動力學(xué)計算過程中,式(6)中質(zhì)量矩陣、外力向量、內(nèi)力向量均進(jìn)行修正,且任意迭代步均可以考慮接觸非線性影響,進(jìn)入接觸范圍的節(jié)點(diǎn)可依據(jù)摩擦因數(shù)進(jìn)行力的傳遞。當(dāng)內(nèi)力I(i)與外力F(i)不平衡時,接觸對將發(fā)生滑動,出現(xiàn)接觸非線性。
采用瞬態(tài)動力學(xué)方法,應(yīng)用商用有限元軟件ABAQUS將擬合的時域試驗條件作為動力學(xué)仿真輸入條件,對閥門進(jìn)行時域振動仿真。邊界條件與假設(shè)如下。
1)考慮閥門線圈螺桿與螺母螺紋之間為接觸關(guān)系,在模型上引入接觸非線性。
2)建立虛擬梁單元兩端分別與閥體墊圈端面、螺母端面耦合約束。
3)旋緊線圈螺桿時,螺紋發(fā)生接觸,螺桿產(chǎn)生位移壓縮梁單元,虛擬梁內(nèi)力即反映螺紋預(yù)緊力,通過預(yù)緊力試驗與螺母旋松旋緊仿真,對螺紋摩擦因數(shù)進(jìn)行精確標(biāo)定。
電磁閥上提取3個響應(yīng)點(diǎn)時域數(shù)據(jù),位置如圖4所示。
圖4 線圈模型Fig.4 Coil model
圖5為經(jīng)過1 s時域隨機(jī)載荷瞬態(tài)動力學(xué)仿真后得到的預(yù)緊力時域曲線,預(yù)緊力由19 kN減小到3.6 kN。預(yù)緊力的減小使得結(jié)構(gòu)整體的剛性減弱,整體模態(tài)頻率降低。
圖5 預(yù)緊力時域曲線Fig.5 Pre-tightening time domain curve
瞬態(tài)動力學(xué)仿真結(jié)果包含所有節(jié)點(diǎn)的時域曲線,如圖6所示,提取線圈根部、中部、尾部3個響應(yīng)點(diǎn)時域響應(yīng)曲線,可見從線圈根部到線圈尾部,加速度響應(yīng)呈現(xiàn)放大趨勢。
圖6 測點(diǎn)時域加速度響應(yīng)Fig.6 Time domain acceleration response of measured point
基于對比分析的需求,應(yīng)用商業(yè)有限元軟件ABAQUS同時進(jìn)行線性模態(tài)仿真。線性模態(tài)仿真基于特征值法,結(jié)構(gòu)自由振動方程為
(7)
式中:C為阻尼矩陣;K為剛度矩陣。
自由振動位移u的振動屬性為
u=φejωt
(8)
模態(tài)計算過程中忽略阻尼C,并考慮式(8),則式(7)轉(zhuǎn)化為
(-ω2M+K)φ=0
(9)
式中ω為結(jié)構(gòu)圓頻率。求解模態(tài)的過程即求解廣義特征值與特征向量。
數(shù)值方法求解特征值過程中,式(9)中無時間相關(guān)量,即質(zhì)量矩陣與剛度矩陣為確定狀態(tài)同時求解多個特征值,從原理上無法考慮與時間相關(guān)的非線性影響,即特征值法無法進(jìn)行接觸非線性模態(tài)分析。
對閥門線性模態(tài)仿真過程中,采用如下邊界條件。
1)假定螺紋為完全約束狀態(tài),即不存在松動、接觸。
2)結(jié)構(gòu)屬性不變,即結(jié)構(gòu)不出現(xiàn)屈服等材料非線性性質(zhì)。
圖7為應(yīng)用線性法得到模態(tài)振型,振型描述了特征值對應(yīng)的特征向量。
圖7 閥門線圈模態(tài)振型Fig.7 Modal vibration of valve coil
模態(tài)頻率對比如表3所示。非線性時域仿真反映出螺紋逐漸松動過程中,預(yù)緊力逐漸變小,結(jié)構(gòu)整體剛性變差,導(dǎo)致模態(tài)頻率降低的結(jié)果,模態(tài)仿真頻率為1.464 kHz;線性隨機(jī)仿真進(jìn)行保守假定,即認(rèn)為螺紋為完全固定狀態(tài),反映在仿真結(jié)果上模態(tài)頻率更高,模態(tài)仿真頻率為1.576 kHz。如圖1所示,推進(jìn)系統(tǒng)振動試驗過程中,電磁閥線圈出現(xiàn)松動,線圈頭部加速度峰值頻率約為1.450 kHz。
表3 模態(tài)頻率對比Tab.3 Comparison of modal frequency
將非線性仿真方法與線性仿真方法結(jié)果進(jìn)行對比分析。將瞬態(tài)動力學(xué)線圈根部、中部、尾部響應(yīng)曲線依據(jù)周期圖法轉(zhuǎn)換為PSD曲線,如圖8(a)所示。圖8(b)為使用模態(tài)疊加法獲得的線圈相應(yīng)位置PSD曲線。將曲線對比可以發(fā)現(xiàn),線性與非線性PSD曲線形態(tài)相似,均反映主振模態(tài)特征,區(qū)別在于瞬態(tài)動力學(xué)采用的時域方法考慮了螺紋松動接觸非線性影響,計算頻率低于線性模態(tài)仿真方法。
圖8 PSD仿真結(jié)果對比Fig.8 Comparison of PSD simulation results
值得說明的是,本文提出的非線性模態(tài)仿真方法無法獲取系統(tǒng)模態(tài)振型,但由于線性、非線性方法獲得的頻率峰值接近,工程仿真過程中將線性模態(tài)振型作為參考。
圖8(a)的非線性仿真方法模型中使用預(yù)緊力試驗標(biāo)定過的摩擦因數(shù),計算過程中考慮了接觸非線性,最終得到模態(tài)非線性仿真結(jié)果;而圖8(b)線性仿真方法無法考慮螺紋副接觸狀態(tài),采用保守的全固定約束方式,得到保守的線性仿真結(jié)果。
雖然線性算法可假定部分螺紋固定,將模態(tài)仿真結(jié)果標(biāo)定至試驗值,但全固定至部分固定變化規(guī)律缺乏理論與物理事實依據(jù),僅可做經(jīng)驗性假設(shè),而這一關(guān)鍵性假設(shè)對仿真結(jié)果影響較大,除完全固定與完全不固定兩種極限狀態(tài)外,無法證明其余狀態(tài)下仿真結(jié)果偏保守或者激進(jìn)。相較而言,在本文采用的非線性仿真方法中,對仿真結(jié)果存在影響的摩擦因數(shù)具備明確物理含義,可采用預(yù)緊力試驗進(jìn)行標(biāo)定,使用標(biāo)定參數(shù)進(jìn)行仿真,得到更為接近試驗值的仿真結(jié)果。
本文采用的非線性模態(tài)仿真方法理論上可推廣應(yīng)用于包括材料非線性、幾何非線性、邊界非線性等在內(nèi)的非線性模態(tài)分析。
本文研究隨機(jī)功率譜的時域數(shù)據(jù)擬合與驗證方法,并應(yīng)用于某閥門隨機(jī)振動線圈松動仿真,得到了閥門非線性模態(tài)頻率,從機(jī)理上定量解釋線圈松動對閥門模態(tài)的非線性影響。
本文提出的非線性模態(tài)仿真方法引入了接觸非線性的影響,通過預(yù)緊力試驗對摩擦因數(shù)進(jìn)行標(biāo)定后,對比采用保守邊界條件的線性方法,仿真結(jié)果與試驗結(jié)果誤差從8.69%降低為0.96%。