考慮彈性影響的乘波體飛行動(dòng)力學(xué)特性

向錦武 曾開春 聶 璐

(北京航空航天大學(xué) 航空科學(xué)與工程學(xué)院，北京100191)

乘波體構(gòu)型近空間高超聲速飛行器(乘波體飛行器)可用于天地往返運(yùn)輸系統(tǒng)和高速導(dǎo)彈，在軍事和民用領(lǐng)域都具有廣闊的應(yīng)用前景.與常規(guī)飛行器相比，乘波體飛行器采用機(jī)體/發(fā)動(dòng)機(jī)一體化設(shè)計(jì)，造成氣動(dòng)/推進(jìn)的強(qiáng)烈耦合［1］.同時(shí)，由于其細(xì)長(zhǎng)體外形的機(jī)體結(jié)構(gòu)中使用了大量輕質(zhì)材料，機(jī)體結(jié)構(gòu)柔性大、振動(dòng)頻率低，在氣動(dòng)載荷激勵(lì)下，可能發(fā)生顯著的彈性變形和彈性振動(dòng).機(jī)體結(jié)構(gòu)的大幅度彈性變形不僅改變了全機(jī)的氣動(dòng)特性，還引起推進(jìn)系統(tǒng)空氣流量以及入流參數(shù)的變化，形成乘波體飛行器特有的氣動(dòng)/彈性/推進(jìn)耦合現(xiàn)象［2－3］.而機(jī)體的低頻彈性振動(dòng)則可能激發(fā)剛體動(dòng)力學(xué)短周期模態(tài)，形成剛?cè)狁詈犀F(xiàn)象.因此，結(jié)構(gòu)彈性加劇了乘波體飛行動(dòng)力學(xué)的耦合性，并使其呈現(xiàn)出多尺度、強(qiáng)非線性的特點(diǎn)，進(jìn)而對(duì)總體設(shè)計(jì)和控制系統(tǒng)設(shè)計(jì)提出了嚴(yán)峻挑戰(zhàn).由于結(jié)構(gòu)彈性參數(shù)會(huì)因設(shè)計(jì)更改而改變，并受到氣動(dòng)加熱等因素的影響，其值可能在較大范圍變化.因此，有必要研究結(jié)構(gòu)彈性參數(shù)對(duì)乘波體飛行器多耦合飛行動(dòng)力學(xué)特性的影響規(guī)律.

乘波體飛行器動(dòng)力學(xué)特性分析的關(guān)鍵是建立多耦合動(dòng)力學(xué)模型.目前，主要的建模方法包括基于FEM(Finite Elenent Melthod)結(jié)構(gòu)模態(tài)計(jì)算的動(dòng)力學(xué)建模［4］、基于 CFD(Computationd Fluid Dynamics)數(shù)值擬合的動(dòng)力學(xué)建模［5］、基于近似理論方法的耦合動(dòng)力學(xué)建模等［6－7］.前2種方法主要面向工程，其分析過程需要詳細(xì)的結(jié)構(gòu)布置與外形尺寸參數(shù)，而近似理論建模方法主要面向科學(xué)研究，其計(jì)算效率相對(duì)較高，是目前乘波體飛行器飛行動(dòng)力學(xué)研究的主要建模方法.因此，本文采用近似理論建模方法，將氣動(dòng)外形、推進(jìn)系統(tǒng)、結(jié)構(gòu)模型等進(jìn)行了簡(jiǎn)化處理，再分別利用近似理論氣動(dòng)力估算方法、瑞利流理論、假設(shè)模態(tài)法等建立了考慮多種耦合因素的動(dòng)力學(xué)模型.研究了在不同機(jī)體結(jié)構(gòu)剛度下，結(jié)構(gòu)彈性對(duì)飛行動(dòng)力學(xué)靜配平特性及動(dòng)態(tài)特性的影響規(guī)律，為總體設(shè)計(jì)及飛行控制系統(tǒng)設(shè)計(jì)提供了參考.

1 高超聲速飛行器飛行動(dòng)力學(xué)模型

典型的乘波體飛行器采用機(jī)體/發(fā)動(dòng)機(jī)一體化設(shè)計(jì).機(jī)體前部呈尖楔形，高超聲速飛行時(shí)，產(chǎn)生附體斜激波，為飛行器提供升力，同時(shí)為發(fā)動(dòng)機(jī)提供高壓空氣;機(jī)體后部是發(fā)動(dòng)機(jī)的外噴管，產(chǎn)生額外的推力，同時(shí)也產(chǎn)生附加的升力和俯仰力矩.因此，動(dòng)力學(xué)建模需要反映氣動(dòng)與推進(jìn)的這種耦合關(guān)系.本文分別采用激波膨脹波方法、瑞利流理論、假設(shè)模態(tài)法等對(duì)氣動(dòng)、推進(jìn)系統(tǒng)和結(jié)構(gòu)彈性動(dòng)力學(xué)進(jìn)行建模，給出了考慮氣動(dòng)/推進(jìn)/結(jié)構(gòu)耦合的乘波體飛行器動(dòng)力學(xué)模型.

1.1 氣動(dòng)力模型

采用基于激波膨脹波方法［7］的當(dāng)?shù)亓骰钊碚摚?］給出高超聲速非定常氣動(dòng)參數(shù)，進(jìn)而建立了乘波體飛行器的氣動(dòng)模型.激波膨脹波方法適用于附體激波流動(dòng)，能給出激波位置，適用的馬赫數(shù)范圍大，是計(jì)算高超聲速附體激波流動(dòng)氣動(dòng)力的基本方法.而牛頓沖擊流理論［6］形式簡(jiǎn)單，可用于快速估算操縱面的氣動(dòng)力，一般適用于馬赫數(shù)5以上.當(dāng)?shù)亓骰钊碚撃軌蚩紤]非定常效應(yīng)，其應(yīng)用范圍為馬赫數(shù)1.2～10，但由于高超聲速流動(dòng)時(shí)的非定常效應(yīng)影響較小，其使用限制可適當(dāng)放寬.當(dāng)?shù)亓骰钊碚撊缡?1)、式(2)所示.

其中，Vn為氣動(dòng)微元瞬時(shí)法向擾動(dòng)速度;q為俯仰角速度;?wd/?x 和分別為彈性變形轉(zhuǎn)角和速度;aSE，pSE，VLocal分別為激波膨脹波方法計(jì)算的當(dāng)?shù)匾羲佟簭?qiáng)和氣流速度模量;利用當(dāng)?shù)亓骰钊碚撔拚膲簭?qiáng)分布，可以積分得到各個(gè)氣動(dòng)面的氣動(dòng)力和力矩.此外，由于考慮了彈性變形的影響，以上建立的氣動(dòng)力模型能夠反映氣動(dòng)與結(jié)構(gòu)彈性的耦合.

1.2 推進(jìn)系統(tǒng)模型

根據(jù)文獻(xiàn)［6］的超燃沖壓發(fā)動(dòng)機(jī)模型，可以將發(fā)動(dòng)機(jī)分成擴(kuò)散段、燃燒室和內(nèi)噴管.在激波膨脹波方法給出上游流動(dòng)參數(shù)和激波位置的基礎(chǔ)上，可以計(jì)算出推進(jìn)系統(tǒng)的入流參數(shù).再利用瑞利增能流［9］理論可以給出內(nèi)流參數(shù)，具體計(jì)算方法參看文獻(xiàn)［6］.在分析內(nèi)流的基礎(chǔ)上，根據(jù)發(fā)動(dòng)機(jī)入流和出口參數(shù)，可給出推力的公式如下:

其中，下標(biāo)i表示入口流動(dòng)參數(shù);下標(biāo)e表示出口參數(shù);下標(biāo)∞表示自由流參數(shù);為發(fā)動(dòng)機(jī)氣體流量.可以發(fā)現(xiàn)，發(fā)動(dòng)機(jī)的狀態(tài)強(qiáng)烈依賴于入流參數(shù)，而入流參數(shù)又決定于上游機(jī)體前部下表面的流動(dòng)參數(shù)與激波位置，因此該建模方法能夠考慮氣動(dòng)與推進(jìn)的耦合.

1.3 機(jī)體結(jié)構(gòu)彈性模型

乘波體飛行器剛度分布中部最大，并向兩端遞減，因此，彈性變形主要集中在兩端.在缺乏數(shù)據(jù)的設(shè)計(jì)初期，這種彈性變形方式適合采用基于彈性集中假設(shè)的兩段模型［6］.該模型在機(jī)身前部和后部連接的中部引入彈性鉸，并假設(shè)變形時(shí)，前部和后部均剛性轉(zhuǎn)動(dòng)，彈性鉸發(fā)生變形.此時(shí)，機(jī)體彈性變形位移wd可表示為模態(tài)廣義坐標(biāo)與模態(tài)函數(shù)的乘積:

其中，由于僅考慮一階模態(tài)，因此i=1.根據(jù)拉格朗日方程可以得出質(zhì)量歸一化的彈性模態(tài)函數(shù)為

其中，J1，J2分別為重心前后機(jī)身段的俯仰轉(zhuǎn)動(dòng)慣量;L1，L2分別為飛行器前后機(jī)身段的長(zhǎng)度.在模態(tài)坐標(biāo)下，引入彈性模態(tài)廣義力來計(jì)算結(jié)構(gòu)彈性變形，其表達(dá)式為

其中，fz為機(jī)體法向分布載荷;Fzj為法向集中載荷.機(jī)體前后部分中線的等效傾角改變量Δ1，Δ2與模態(tài)廣義坐標(biāo)的關(guān)系為

高超聲速飛行器巡航飛行時(shí)，氣動(dòng)加熱對(duì)結(jié)構(gòu)力學(xué)特性影響較大，但相對(duì)于氣動(dòng)彈性響應(yīng)，這種影響是一個(gè)緩慢的過程，因此本文將氣動(dòng)加熱的影響考慮成結(jié)構(gòu)剛度在一定范圍內(nèi)變化.

2 縱向飛行力學(xué)特性分析方法

從基本的拉格朗日方程出發(fā)，推導(dǎo)了縱向剛?cè)狁詈系娘w行動(dòng)力學(xué)模型.機(jī)體坐標(biāo)系下，柔性飛行器縱向運(yùn)動(dòng)學(xué)和耦合動(dòng)力學(xué)方程如下:

其中，狀態(tài)變量h為高度;u為x向飛行速度;w為z向速度;θ為俯仰角;q為俯仰角速度;η為彈性模態(tài)廣義坐標(biāo)為彈性模態(tài)廣義速度;下文分析中使用的迎角α可以近似由w與u的比值給出;Fx，F(xiàn)z分別為x向和z向合力;M為俯仰合力矩;N為彈性廣義力.

2.1 靜態(tài)配平與飛行包線

求解靜配平問題時(shí)，需要首先對(duì)方程式(8)～式(13)進(jìn)行處理，即令飛行速度分量的時(shí)間導(dǎo)數(shù)和、俯仰角速度 q及其對(duì)時(shí)間的導(dǎo)數(shù)、彈性變形速度均為0，將原來的6個(gè)方程化簡(jiǎn)為5個(gè).然后在指定飛行高度和速度，利用Newton-Raphson方法求解關(guān)于迎角α、俯仰角θ、升降舵偏轉(zhuǎn)角δe(向下偏轉(zhuǎn)為正)、發(fā)動(dòng)機(jī)油門δT以及彈性廣義坐標(biāo)η等5個(gè)變量的非線性方程組，可獲得靜平衡解.

對(duì)于飛行包線，考慮以下3種邊界條件:

1)最大動(dòng)壓邊界:熱防護(hù)系統(tǒng)的耐受范圍決定了飛行的最大動(dòng)壓，其準(zhǔn)確值需要由試驗(yàn)才能給出，本文取為95.76 kPa(2000 lb/ft).

2)發(fā)動(dòng)機(jī)熱阻塞邊界:當(dāng)發(fā)動(dòng)機(jī)燃燒室出口馬赫數(shù)低于1時(shí)，發(fā)動(dòng)機(jī)會(huì)因熱阻塞而不能正常工作.該邊界是由可壓縮流體力學(xué)理論［9］推導(dǎo)得出的，本文給出一定的裕度，即要求發(fā)動(dòng)機(jī)燃燒室出口馬赫數(shù)不得低于1.3.

3)最大油門邊界:由于進(jìn)入發(fā)動(dòng)機(jī)的空氣有限，因此當(dāng)其中的氧氣完全與燃油發(fā)生化學(xué)反應(yīng)，即等效油門為1后，不能通過增大油門來增加推力.同時(shí)考慮到進(jìn)入發(fā)動(dòng)機(jī)的氧氣并不能都發(fā)生化學(xué)反應(yīng)，因此等效油門最大值取為0.9.

上述的3個(gè)邊界條件中，條件1)確定了可配平飛行的最大速度和最低高度，其與飛行器參數(shù)無關(guān)，僅取決于大氣環(huán)境和飛行速度;條件2)和條件3)確定了最低飛行速度和最大高度，二者均與飛行器的參數(shù)相關(guān).本文在求解條件2)和條件3)時(shí)，首先選取了用于分析求解的馬赫數(shù)值;然后分別針對(duì)每個(gè)馬赫數(shù)值，使用二分法求解最大可配平高度;最后在速度-高度圖中，將各個(gè)離散的邊界點(diǎn)連接起來形成飛行包線.

2.2 動(dòng)態(tài)特性

在配平狀態(tài)下對(duì)方程式(8)～式(13)做線性化處理，可以得到如式(14)、式(15)所示的線性動(dòng)力學(xué)分析模型.

3 機(jī)體剛度對(duì)飛行力學(xué)特性的影響

根據(jù)文獻(xiàn)［7］給出的典型乘波體飛行器參數(shù)，從飛行動(dòng)力學(xué)靜配平特性、飛行包線要求、動(dòng)態(tài)特性3個(gè)角度研究了機(jī)體結(jié)構(gòu)剛度對(duì)飛行動(dòng)力學(xué)特性的影響規(guī)律.

3.1 機(jī)體結(jié)構(gòu)剛度對(duì)靜態(tài)特性的影響

機(jī)體一階彎曲自然頻率ω1代表了機(jī)體彎曲剛度的大小，其值越大表明機(jī)體彎曲剛度越大.以下分析了配平參數(shù)隨ω1的變化關(guān)系，其中ω1為圓頻率，單位為rad/s.

圖1 前體彈性變形等效傾角隨馬赫數(shù)的變化關(guān)系

圖2 迎角配平值隨馬赫數(shù)的變化關(guān)系

圖3 配平狀態(tài)下的全機(jī)升阻比隨馬赫數(shù)的變化關(guān)系

圖4 油門配平值比隨馬赫數(shù)的變化關(guān)系

圖5 升降舵配平值比隨馬赫數(shù)的變化關(guān)系

圖1～圖5為H=25 908 m(85 000 ft)，Ma=8，ω1分別為無窮大(剛體)，25，18，14 rad/s時(shí)，定常平飛狀態(tài)參數(shù)和控制變量隨馬赫數(shù)的變化.可以看出:

1)在高馬赫數(shù)下，結(jié)構(gòu)彈性對(duì)飛行器的配平狀態(tài)有較大的影響.對(duì)于同一機(jī)體剛度條件下的柔性飛行器，隨著馬赫數(shù)的增大，結(jié)構(gòu)彈性變形引起的機(jī)體前部?jī)A角急劇增大，甚至與設(shè)計(jì)角度接近，因而對(duì)靜配平特性產(chǎn)生的影響不可忽略.

2)機(jī)體剛度下降會(huì)導(dǎo)致配平迎角減小.從圖1、圖2可以看出，在高馬赫數(shù)段，隨著結(jié)構(gòu)剛度減小，結(jié)構(gòu)變形引起機(jī)身前部下表面傾角增大，進(jìn)而導(dǎo)致所需的配平迎角急劇減小.

3)機(jī)體剛度下降會(huì)導(dǎo)致飛行性能降低.從圖3、圖4可以看出，隨著機(jī)體剛度下降，飛行器定常飛行時(shí)的升阻比將減小，所需的油門增大，并且在高馬赫數(shù)段，該影響變得更加顯著.

4)機(jī)體剛度下降會(huì)導(dǎo)致升降舵配平規(guī)律改變.從圖5可以發(fā)現(xiàn)，在所有馬赫數(shù)段，剛體飛行器的升降舵配平值隨馬赫數(shù)增大而減小，為正操縱;而對(duì)于ω1＜25 rad/s的柔性飛行器，在低馬赫數(shù)區(qū)段，升降舵配平值隨馬赫數(shù)增大而增大，為正操縱;在高馬赫數(shù)區(qū)段，升降舵配平值隨馬赫數(shù)增大而減小，變?yōu)榉床倏v.因此控制系統(tǒng)設(shè)計(jì)時(shí)，有必要檢驗(yàn)柔性高超聲速飛行器是否存在高馬赫數(shù)區(qū)反操縱現(xiàn)象.

圖6 是 ω1分別為∞，18，14，10rad/s時(shí)，乘波體飛行器的飛行包線.可以看出機(jī)體剛度減小會(huì)顯著縮小高馬赫數(shù)段的可配平區(qū)域.隨著ω1的減小，發(fā)動(dòng)機(jī)熱阻塞邊界向馬赫數(shù)升高方向移動(dòng)，而最大油門邊界向高度減小方向移動(dòng)，隨著馬赫數(shù)升高，其減小量增大.此外，還可以發(fā)現(xiàn)在 ω1=10 rad/s時(shí)，在馬赫數(shù)16附近，最大油門邊界存在不連續(xù)現(xiàn)象.分析其產(chǎn)生的機(jī)理首先要從高度對(duì)配平特性的影響入手，限于篇幅，本文直接給出結(jié)論:高度降低，剛體飛行器配平迎角減小，配平阻力減小，并且結(jié)構(gòu)彈性變形增大.在該結(jié)論的基礎(chǔ)上，可以發(fā)現(xiàn)，由于結(jié)構(gòu)變形引起阻力增大，因此，當(dāng)彈性變形引起的阻力增加與配平迎角減小引起的阻力減小相抵消時(shí)，通過降低高度的方式不能配平飛行器.

圖6 不同機(jī)身彎曲剛度下的飛行包線

3.2 機(jī)體剛度對(duì)動(dòng)態(tài)特性的影響

在H=25908 m，Ma=8狀態(tài)下，分別計(jì)算ω1為∞，25，18，14 rad/s，升降舵為輸入，軌跡角為輸出的系統(tǒng)開環(huán)的零極點(diǎn)，如圖7所示.可以看出，高超聲速飛行器縱向飛行動(dòng)力學(xué)存在3個(gè)模態(tài)，由5個(gè)特征根組成.其中短周期特征根由絕對(duì)值較大的一正實(shí)根和一負(fù)實(shí)根組成，表明其短周期缺少回復(fù)力，為靜發(fā)散.隨著機(jī)體彎曲剛度的減小，正實(shí)根進(jìn)一步增大，加劇了短周期的不穩(wěn)定.長(zhǎng)周期模態(tài)由2個(gè)小阻尼、低頻率共軛復(fù)根組成，隨著機(jī)體彎曲剛度減小，其實(shí)部變化較小，而虛部向?qū)嵼S靠近.高度模態(tài)對(duì)應(yīng)一個(gè)絕對(duì)值接近0的實(shí)根，其隨機(jī)體剛度變化較小.此外，還可以看出，機(jī)體彎曲剛度對(duì)零點(diǎn)位置的影響較小.綜上所述，機(jī)體剛度降低，短周期模態(tài)穩(wěn)定性急劇下降，進(jìn)而使按照剛體模型設(shè)計(jì)的控制系統(tǒng)失效.

圖7 開環(huán)零極點(diǎn)隨機(jī)身剛度減小的變化

4 結(jié)論

通過對(duì)考慮機(jī)體彈性的乘波體飛行器的縱向飛行動(dòng)力學(xué)建模、計(jì)算，分析了靜態(tài)配平特性、飛行包線及動(dòng)穩(wěn)定性隨機(jī)體彎曲剛度的變化規(guī)律，結(jié)果表明:

1)隨著機(jī)體剛度減小，配平升阻比降低，耗油率增大，由發(fā)動(dòng)機(jī)熱阻塞和最大油門確定的飛行包線縮小;

2)當(dāng)機(jī)體剛度較小時(shí)，飛行器升降舵配平規(guī)律在不同馬赫數(shù)區(qū)段將發(fā)生變化，在低馬赫數(shù)段升降舵操縱為正操縱，在高馬赫數(shù)段變?yōu)榉床倏v;

3)縱向存在不穩(wěn)定的短周期模態(tài)，機(jī)體剛度減小，長(zhǎng)周期模態(tài)頻率將降低，短周期模態(tài)正實(shí)特征根增大，從而使系統(tǒng)穩(wěn)定性降低.

References)

［1］Fidan B，Mirmirani M，Ioannou P.Flight dynamics and control of air-breathing hypersonic vehicles:review and new Directions［R］.AIAA-2003-708，2003

［2］Cockrell C E，Engelund W C.Integrated aeropropulsive computational fluid dynamics methodology for the hyper-X flight experiment［J］.Journal of Spacecraft and Rockets，2001，38(6):836 －843

［3］黃喜元，王青，董朝陽.基于動(dòng)態(tài)逆的高超聲速飛行器魯棒自適應(yīng)控制［J］.北京航空航天大學(xué)學(xué)報(bào)，2011，37(5):560－563 Huang Xiyuan，Wang Qing，Dong Chaoyang.Robust adaptive control for hypersonic vehicle based on dynamic inversion［J］.Journal of Beijing University of Aeronautics and Astronautics，2011，37(5):560 －563(in Chinese)

［4］Heeg J，Zeiler T A，Pototzky AS，et al.Aerothermoelastic analysis of a NASP demonstrator model［R］.NASA-TM-109007，1993

［5］Clark A D，Mirmirani M D，Wu C，et al.An aero-propulsion integrated elastic model of a generic airbreathing hypersonic vehicle［R］.AIAA-2006-6560，2006

［6］Chavez F R，Schmidt D K.Analytical Aeropropulsive/Aeroelastic hypersonic vehicle model with dynamic analysis［J］.AIAA Journal of Guidance，Control，and Dynamics，1994，17(6):1308 －1319

［7］Bolender M A，Doman D B.Nonlinear longitudinal dynamical model of an air-Breathing hypersonic vehicle［J］.Journal of Spacecraft and Rockets，2007，44(2):374 －387

［8］McNamara J J，Crowell A R，F(xiàn)riedmann P R，et al.Approximate modeling of unsteady aerodynamics for hypersonic aeroelasticity［J］.Joural of Aircraft，2010，47(6):1932 －1945

［9］Anderson J D.Modern compressible flow［M］.2nd ed.New York:McGraw-Hill College，1995:77 －91

［10］方振平，陳萬春，張曙光.航空飛行器飛行動(dòng)力學(xué)［M］.北京:北京航空航天大學(xué)出版社，2005:219－220 Fang Zhenping，Chen Wanchun，Zhang Shuguang.Aircraft flight dynamics［M］.Beijing:Beihang University Press，2005:219 －220(in Chinese)