崔建峰, 張科, 呂梅柏
(1.西北工業(yè)大學(xué) 航天學(xué)院, 陜西 西安 710072; 2.航天飛行動力學(xué)國家級重點(diǎn)實(shí)驗(yàn)室, 陜西 西安 710072)
高超聲速飛行器與傳統(tǒng)飛行器相比,在結(jié)構(gòu)上普遍采用輕質(zhì)柔性材料,其氣動外形一般為細(xì)長體、升力體布局、完全或部分乘波體布局。特殊的結(jié)構(gòu)選擇和氣動布局使得飛行器結(jié)構(gòu)的固有振動頻率降低,剛體模態(tài)與彈性模態(tài)的耦合問題更為突出[1],這需要在高超聲速飛行器建模時(shí)考慮飛行器模型為彈性體。
目前,在彈性高超聲速飛行器控制系統(tǒng)研究中,得到廣泛應(yīng)用的是由Bolender、Doman等人提出的二維高超聲速飛行器彈性體動力學(xué)模型[2-3]。國內(nèi)一些學(xué)者也借鑒此種解析建模方法構(gòu)建了相應(yīng)的吸氣式高超聲速飛行器二維動力學(xué)模型[4-5]。上述模型在建模過程中對非定常氣動力的計(jì)算一般以機(jī)身前、后部彈性變形頂端偏轉(zhuǎn)角Δτi引起的氣動力變化來估算非定常氣動力,使用等截面梁估算飛行器的彈性模態(tài)。這種估算方法對于一階振型對應(yīng)的非定常氣動力有較好的近似,但對于更高階振型對應(yīng)的非定常氣動力計(jì)算將產(chǎn)生較大的誤差。另外,經(jīng)后續(xù)研究發(fā)現(xiàn),變截面慣性矩對飛行器的振型模態(tài)影響較大[6],這對飛行器的非定常氣動力計(jì)算會產(chǎn)生較大的影響。
上述模型所提供的氣動數(shù)據(jù)離散,不便于進(jìn)行非線性控制系統(tǒng)綜合設(shè)計(jì)。因此,Parker等人在解析模型的基礎(chǔ)上對氣動數(shù)據(jù)進(jìn)行曲線擬合,得到面向控制的曲線擬合模型[7-8](CFM:curve-fitted model)。氣動數(shù)據(jù)曲線擬合是一個(gè)較為耗費(fèi)時(shí)間與計(jì)算機(jī)資源的過程。所以擬合前,需要依據(jù)可供使用的計(jì)算機(jī)資源選取適當(dāng)簡化形式的模型,這將限制擬合模型的可信度及其與原始物理模型之間的一致程度。因此,建立一種快速、可靠的擬合建模方法十分必要。
針對上述問題,本文在現(xiàn)有彈性模型基礎(chǔ)上,使用有限元法基于變截面自由梁計(jì)算高超聲速飛行器的結(jié)構(gòu)彈性模型,利用當(dāng)?shù)亓骰钊碚撚?jì)算彈性變形引起的非定常氣動力;然后借助均勻設(shè)計(jì)、逐步回歸等統(tǒng)計(jì)學(xué)方法快速高效地獲取相應(yīng)的CFM模型;最后仿真驗(yàn)證模型的可靠性。
本文以文獻(xiàn)[2-3]中給出的高超聲速飛行器模型為主要參考研究對象。該飛行器幾何結(jié)構(gòu)參見圖1,外形參數(shù)詳細(xì)信息可查看文獻(xiàn)[2-3]。飛行器質(zhì)量為11×103kg/m,轉(zhuǎn)動慣量4.45×105kg·m2/m??紤]加熱后材料結(jié)構(gòu)平均彈性模量為7.33×107Pa。
圖1 高超聲速飛行器模型幾何結(jié)構(gòu)
飛行器定常氣動力與推力計(jì)算方法簡述如下:利用空氣動力學(xué)中的斜基波-膨脹波理論分別在高超聲速飛行器前體斜坡、發(fā)動機(jī)燃燒室、后體膨脹面以及控制舵面處計(jì)算飛行器表面壓強(qiáng)分布及氣流參數(shù);由于波后壓強(qiáng)分布均勻,從而可得到各面元上的定常氣動力和力矩;飛行器超燃沖壓發(fā)動機(jī)的推力則通過瑞利流模型進(jìn)行估算[2-3]。
對于高超聲速飛行器彈性模型而言,更為關(guān)鍵的是非定常氣動力的有效計(jì)算,因此,本文主要對飛行器的彈性模態(tài)、非定常氣動力求解做詳細(xì)介紹。
對于變截面非均勻歐拉-貝努利,其無阻尼自由振動方程可表示為[9]:
(1)
為計(jì)算變截面自由梁的振型模態(tài)及相應(yīng)的模態(tài)頻率,在文獻(xiàn)[6]中,采用了“假設(shè)模態(tài)法”。假設(shè)模態(tài)法可以較為精確地求解模態(tài)頻率,但由于初始假設(shè)模態(tài)的偏差,通過一次計(jì)算難以保證模態(tài)特征向量的估算精度。因此,如需獲取更為精確的振型,可以使用子空間迭代法或有限元方法進(jìn)行計(jì)算[9]。其中,有限元法可以看成是分區(qū)的瑞利-里茲法。其首先將一個(gè)復(fù)雜結(jié)構(gòu)分成多個(gè)區(qū)域,然后在每個(gè)區(qū)域內(nèi)使用瑞利-里茲法來處理,最后通過邊界條件聯(lián)合求解。經(jīng)過這樣的剖分后,各單元幾何形狀簡單規(guī)則,單元的初始假設(shè)模態(tài)易于選取,且偏差較小,可以明顯提高最終求解精度。本文將使用有限元方法求解飛行器等效變截面自由梁的振型模態(tài)。
由于參考研究對象為二維飛行器,因此,本文在此進(jìn)行如下假設(shè):二維飛行器可擴(kuò)展為單位寬度的三維飛行器,三維擴(kuò)展飛行器橫截面的截面形式為長方形。因此,截面慣性矩根據(jù)公式(2)進(jìn)行計(jì)算。
Iy=bh3/12+bhz2
(2)
式中:b為截面寬度,h為截面高度;z為截面形心主軸與機(jī)體坐標(biāo)系x軸之間的距離。
在有限元分析中,選擇梁單元的類型為三維有限應(yīng)變梁單元BEAM189,設(shè)置截面類型為RECT,并根據(jù)前述參數(shù)設(shè)置材料屬性。在定義位移邊界條件后,使用分塊法模態(tài)分析求解器求解。經(jīng)過有限元計(jì)算獲取飛行器前三階模態(tài)頻率為:15.5 rad/s,29.1 rad/s,44.8 rad/s,前三階模態(tài)振型如圖2所示。
圖2 質(zhì)量歸一化彈性模態(tài)振型
在圖2中,同時(shí)給出了等截面自由等效梁、等截面懸臂等效梁及使用假設(shè)模態(tài)法計(jì)算的變截面梁的一階振型模態(tài)。從圖中對比可以看出,變截面慣性矩對飛行器的振型模態(tài)影響較大,尤其在飛行器頭部0~5 m位置,明顯大于等截面梁的振型變化;但在計(jì)算對象同為變截面梁的情況下,假設(shè)模態(tài)法對振型的估算與更為精確的有限元方法相比,仍存在較大的差距。如果使用等截面梁的振型模態(tài)或假設(shè)模態(tài)法估算的振型模態(tài)來計(jì)算飛行器的非定常氣動力必然會產(chǎn)生較大誤差。另外,通過圖中的變截面梁高階模態(tài)振型曲線可以直觀看出,如果使用等效攻角增量來估算非定常氣動力,由于高階模態(tài)振型中部隆起,必然導(dǎo)致計(jì)算結(jié)果存在偏差。
在獲取高超聲速飛行器的彈性模態(tài)后,本文使用斜基波-膨脹波理論所計(jì)算的定常流場,基于當(dāng)?shù)亓骰钊碚撚?jì)算非定常氣動力。
在經(jīng)典活塞理論的基礎(chǔ)上保留一階項(xiàng)并將來流的參數(shù)替換為當(dāng)?shù)亓鲃訁?shù)(用下標(biāo)L表示)后,物面壓強(qiáng)可表示為[10]:
(3)
式中:n0為物面變型前的外法線單位矢量,n為物面變型后的外法線單位矢量,VL和Vb分別表示當(dāng)?shù)亓鲃铀俣群臀锩嬲駝铀俣?W是當(dāng)?shù)叵孪此俣?由物面變形VL·δn和振動Vb·n組成。
運(yùn)用該方法求解非定常氣動力時(shí),基于定常流場得到物面當(dāng)?shù)氐牧鲃訁?shù),再通過(3)式算出物面運(yùn)動時(shí)當(dāng)?shù)叵孪此俣鹊淖兓?從而得到物面壓強(qiáng)隨時(shí)間的變化,將壓強(qiáng)沿物面積分即可得到非定常氣動力。
利用解析方法求解的各狀態(tài)下的氣動力與推力是離散數(shù)據(jù),不便于進(jìn)行非線性控制系統(tǒng)的綜合設(shè)計(jì),有必要對離散氣動數(shù)據(jù)進(jìn)行曲線擬合,得到各力與飛行器狀態(tài)之間的近似函數(shù)關(guān)系。傳統(tǒng)的擬合方法首先需要獲取大量的擬合數(shù)據(jù)集和驗(yàn)證數(shù)據(jù)集驗(yàn)證,這個(gè)過程較為耗費(fèi)時(shí)間與計(jì)算機(jī)資源,且因數(shù)據(jù)集龐大,給后期模型擬合帶來不便。因此,本文使用均勻試驗(yàn)設(shè)計(jì)來簡化數(shù)據(jù)生成的過程。
在確定試驗(yàn)因素及試驗(yàn)設(shè)計(jì)空間后,需要根據(jù)均勻設(shè)計(jì)的原理和方法構(gòu)造指定試驗(yàn)因素的均勻設(shè)計(jì)表。
常用的均勻設(shè)計(jì)表構(gòu)造方法有好格子點(diǎn)法、方冪好格子點(diǎn)法及在它們基礎(chǔ)上進(jìn)行改進(jìn)的修正法與切割法[11]。因本文中的試驗(yàn)因素均為連續(xù)變量,對變量水平?jīng)]有嚴(yán)格限制,故本文使用好格子點(diǎn)法、方冪好格子點(diǎn)法及其修正方法構(gòu)建均勻設(shè)計(jì)表。詳細(xì)的方法與步驟請參見文獻(xiàn)[11]。
表1 狀態(tài)變量工作范圍
均勻設(shè)計(jì)表的構(gòu)造過程中,需要在不同生成表間進(jìn)行尋優(yōu)。而比較2個(gè)均勻設(shè)計(jì)表的好壞就等價(jià)于比較它們所對應(yīng)2組點(diǎn)集的均勻性。常用的均勻性測度有Lp-星偏差、中心化L2-偏差(CD2)等[11]。根據(jù)上述均勻性測度的準(zhǔn)確性及計(jì)算復(fù)雜程度,本文選用中心化L2-偏差作為均勻設(shè)計(jì)表的均勻性測度:
以CD2(P)≤0.12為評價(jià)標(biāo)準(zhǔn),最終獲取的均勻設(shè)計(jì)表為U*(240,24012),其中星號表示由修正方法生成,該表的中心化L2-偏差CD2(P)=0.117 2。
利用上節(jié)獲取的均勻設(shè)計(jì)表U*(240,24012)進(jìn)行仿真試驗(yàn),生成相應(yīng)的氣動力、推力及廣義力數(shù)據(jù)。在這些數(shù)據(jù)的基礎(chǔ)上,為獲取各變量間的相關(guān)關(guān)系,需要進(jìn)行回歸分析。根據(jù)理論及實(shí)際計(jì)算結(jié)果分析可知,各力與試驗(yàn)因素間并非簡單的線性關(guān)系,需要考慮到各因素之間的交互作用及因素高次項(xiàng)的影響。但將如此多的因素全部引入到擬合模型中,勢必過于復(fù)雜,且可能形成“過擬合”。因此,需要從這些眾多相關(guān)因素中挑選出對結(jié)果有顯著影響的部分因素,剔除一些不顯著因子,故可以使用逐步回歸分析法。逐步回歸方法具體步驟可參考文獻(xiàn)[11]。
本文最終獲取的CFM模型為:
(4)
(4)式中,角度單位為rad;空氣密度則根據(jù)美國1976年標(biāo)準(zhǔn)大氣表計(jì)算得到,具體計(jì)算方法可查看相關(guān)文獻(xiàn)。從(4)式中可以看出,飛行器升力及俯仰力矩系數(shù)中彈性廣義坐標(biāo)及其速度相關(guān)項(xiàng)反映的是非定常氣動力系數(shù),其與相應(yīng)的氣動力與力矩系數(shù)(包括廣義力系數(shù))基本呈線性關(guān)系。但相比于剛體參數(shù)項(xiàng),彈性項(xiàng)影響較弱。阻力系數(shù)及推力系數(shù)因受剛體參數(shù)影響更為顯著,因此,在逐步回歸過程中,與之有關(guān)的彈性氣動導(dǎo)數(shù)項(xiàng)作為非顯著因子被剔除。
在獲取CFM模型后,需要對模型進(jìn)行驗(yàn)證。一般可以通過復(fù)相關(guān)系數(shù)R與殘余標(biāo)準(zhǔn)偏差S來檢驗(yàn)回歸方程的顯著性與精度[11]。
對于驗(yàn)證數(shù)據(jù)集的生成,同樣可以使用均勻設(shè)計(jì)表來獲取。在獲取驗(yàn)證數(shù)據(jù)集后,本文CFM模型各統(tǒng)計(jì)量數(shù)值如表2所示:
表2 CFM模型統(tǒng)計(jì)量
從表2中可以看出,飛行器各力回歸方程的復(fù)相關(guān)系數(shù)均接近于1,說明CFM模型與原物理模型一致程度較高。
對于已有的CFM模型,文獻(xiàn)[8]在吸氣式高超聲速飛行器幾何外形與配置上與本文的基本一致。但由于非定常氣動力計(jì)算方法不一致,因此,本文中僅與其比較剛體模型相關(guān)氣動力、推力及俯仰力矩的擬合情況。為了直觀地查看擬合效果,圖3給出了使用真實(shí)物理模型、本文CFM模型、文獻(xiàn)CFM模型計(jì)算得到氣動力與力矩隨速度與攻角的變化圖。
圖3 氣動力與力矩?cái)M合對比
從圖3中可以看出,相比于文獻(xiàn)[8]中的CFM模型,本文借用均勻設(shè)計(jì)法所獲取的CFM模型與真實(shí)物理模型更為一致,模型精度也更好。但用于擬合及驗(yàn)證的數(shù)據(jù)點(diǎn)卻不到2 000個(gè),建模過程耗費(fèi)的時(shí)間與計(jì)算機(jī)資源大為減少。
由于非定常氣動力算法的改進(jìn),系統(tǒng)的平衡和動態(tài)特性將發(fā)生變化。為了研究飛行器平衡和動態(tài)特性,選取1組高超聲速飛行器典型飛行條件:(2 100 m/s,22 km)、(2 100 m/s,32 km)、(3 100 m/s,22 km)(3 100 m/s,32 km),然后分別在原有模型與改進(jìn)模型下,得到從(δe,φ)至(V,h)的飛行器傳輸零極點(diǎn)分布,如圖4所示。
圖4 零極點(diǎn)分布對比圖
通過圖4的對比可以看出,相對原有高超聲速飛行器彈性模型,使用本文方法構(gòu)建的彈性模型動態(tài)特性有較大的改變。原有弱不穩(wěn)定的彈性模態(tài)(3組共軛零極點(diǎn))在計(jì)入變截面影響并改進(jìn)非定常氣動力算法后,將呈現(xiàn)出部分不穩(wěn)定的特征,較為顯著的是一階彈性模態(tài)的變化,代表其特征的極點(diǎn)由s左半平面進(jìn)入右半平面,如圖5中局部放大圖所示。因此在高超聲速飛行器控制系統(tǒng)設(shè)計(jì)中必需要對彈性模態(tài)進(jìn)行穩(wěn)定,否則所設(shè)計(jì)的控制系統(tǒng)必因彈性模態(tài)發(fā)散而失控。
本文提出了一種面向控制的彈性高超聲速飛行器建模方法。該方法對現(xiàn)有模型中結(jié)構(gòu)彈性模型的構(gòu)建及非定常氣動力的計(jì)算方法進(jìn)行了改進(jìn),使用變截面自由梁來構(gòu)建高超聲速飛行器的結(jié)構(gòu)彈性模型,利用當(dāng)?shù)亓骰钊碚撚?jì)算彈性變形引起的非定常氣動力;并借助均勻設(shè)計(jì)、逐步回歸等統(tǒng)計(jì)學(xué)理論構(gòu)建了彈性高超聲速飛行器的CFM模型。通過與現(xiàn)有CFM模型的各項(xiàng)數(shù)據(jù)對比分析表明,本文中的彈性高超聲速飛行器建模方法高效、可靠,其所建立的CFM模型與原始物理模型的一致程度及精度均優(yōu)于現(xiàn)有CFM模型,而所消耗的時(shí)間與計(jì)算機(jī)資源遠(yuǎn)小于現(xiàn)有擬合方法。另外,通過對所構(gòu)建模型進(jìn)行的動態(tài)特性分析,指出了高超聲速飛行器彈性模型將呈現(xiàn)出不穩(wěn)定的彈性模態(tài)特征。這需要在后續(xù)控制器設(shè)計(jì)中引起關(guān)注。
參考文獻(xiàn):
[1] 楊超,許赟,謝長川. 高超聲速飛行器氣動彈性力學(xué)研究綜述[J]. 航空學(xué)報(bào),2010, 31(1): 1-11
Yang Chao, Xu Yun, Xie Changchuan. Review of Studies on Aero-Elasticity of Hypersonic Vehicles[J]. ACTA Aeronautica et Astronautica Sinica, 2010, 31(1): 1-11 (in Chinese)
[2] Bolender M A, Doman D B. A Non-Linear Model for the Longitudinal Dynamics of a Hypersonic Air-Breathing Vehicle[C]∥San Francisco, CA, United States: Collection of Technical Papers-AIAA Guidance, Navigation, and Control Conference, 2005: 3937-3958
[3] Bolender M A, Doman D B. Flight Path Angle Dynamics of Air-Breathing Hypersonic Vehicles[C]∥AIAA Guidance, Navigation, and Control Conference and Exhibit, Keystone, Colorado, 2006
[4] 孟中杰,陳凱,黃攀峰,等. 高超聲速飛行器機(jī)體/發(fā)動機(jī)耦合建模與控制[J]. 宇航學(xué)報(bào). 2008, 29(5): 1509-1514
Meng Zhongjie, Chen Kai, Huang Panfeng, et al. The Coupling Model and Control between Scramjet and Airframe for Hypersonic Vehicle [J]. Journal of Astronautics. 2008, 29(5):1509-1514 (in Chinese)
[5] 張冉,李惠峰,肖進(jìn). 高超聲速飛行器剛體/彈性體耦合動力學(xué)建模[J]. 北京航空航天大學(xué)學(xué)報(bào), 2012(02): 160-165
Zhang Ran, Li Huifeng, Xiao Jin. Hypersonic Vehicle Rigid/Elastic Coupled Dynamic Modeling [J]. Journal of Beijing University of Aeronautics and Astronautics, 2012(02): 160-165 (in Chinese)
[6] 蘇二龍,羅建軍,黃興李,等. 考慮氣動加熱和變截面慣性矩的高超聲速飛行器建模與分析[J]. 宇航學(xué)報(bào), 2012, 33(6): 690-697
Su Erlong, Luo Jianjun, Huang Xingli, et al. Modeling and Analysis of Hypersonic Vehicle Considering Variable Cross-Section Moment of Inertia and Aerodynamics Heating [J]. Journal of Astronautics, 2012, 33(6): 690-697 (in Chinese)
[7] Parker J T, Serrani A, Yurkovich S, et al. Control-Oriented Modeling of an Air-Breathing Hypersonic Vehicle[J]. Journal of Guidance, Control, and Dynamics, 2007, 30(3): 856-869
[8] Fiorentini L. Nonlinear Adaptive Controller Design for Air-Breathing Hypersonic Vehicles [D]. Ohio: The Ohio State University, 2010: 122-159
[9] 諸德超,程偉,邢譽(yù)峰,等. 工程振動基礎(chǔ)[M]. 北京:北京航空航天大學(xué)出版社, 2004
Zhu D C, Cheng W, Xing Yufeng, et al. Foundation of Engineering Vibration[M]. Beijing: Beihang University Press. 2004: 122-159 (in Chinese)
[10] 張偉偉, 史愛明, 王剛,等. 結(jié)合定常CFD 技術(shù)的當(dāng)?shù)亓骰钊碚揫J]. 西北工業(yè)大學(xué)學(xué)報(bào), 2004, 22(5): 545-549
Zhang Weiwei, Shi Aiming, Wang Gang, et al. On Determining Unsteady Aerodynamic Loads Accurately and Efficiently[J]. Journal of Northwestern Polytechnical University, 2004, 22(5): 545-549 (in Chinese)
[11] 方開泰,劉民千,周永道. 試驗(yàn)設(shè)計(jì)與建模[M]. 北京:高等教育出版社, 2011: 148-201
Fang K T, Liu M Q, Zhou Y D. Design and Modeling of Experiments[M]. Beijing: Higher Education Press, 2011: 148-201 (in Chinese)