王 斐，冉玉國(guó)，李秋彥

航空工業(yè)成都飛機(jī)設(shè)計(jì)研究所，四川 成都 610091

引言

氣動(dòng)彈性是一門多學(xué)科交叉的力學(xué)分支[1-2]，在傳統(tǒng)飛機(jī)研制流程中，由于前期飛機(jī)結(jié)構(gòu)布置、尺寸參數(shù)以及重量特性等輸入條件不具備或存在諸多不確定性，氣動(dòng)彈性專業(yè)介入一般相對(duì)滯后。而在飛機(jī)研制后期，一旦暴露氣動(dòng)彈性穩(wěn)定性問(wèn)題，可設(shè)計(jì)空間有限，改進(jìn)優(yōu)化難度大，甚至可能導(dǎo)致設(shè)計(jì)方案被推翻等災(zāi)難性問(wèn)題，付出人力、時(shí)間等巨大的代價(jià)[3]。

現(xiàn)代飛機(jī)由于不斷追求減輕結(jié)構(gòu)重量和提高飛行性能的目標(biāo)，使得氣動(dòng)彈性穩(wěn)定性問(wèn)題變得愈加突出。這就要求氣動(dòng)彈性專業(yè)設(shè)計(jì)人員在飛機(jī)研制方案設(shè)計(jì)階段盡早開(kāi)展工作，避免對(duì)后續(xù)設(shè)計(jì)帶來(lái)顛覆性問(wèn)題。飛機(jī)研制是一項(xiàng)復(fù)雜的系統(tǒng)工程，在方案設(shè)計(jì)階段，總體外形和結(jié)構(gòu)布局不斷迭代更新，沒(méi)有足夠的周期及成本用于對(duì)各個(gè)方案開(kāi)展詳細(xì)設(shè)計(jì)，不可能具備具體結(jié)構(gòu)參數(shù)和重量數(shù)據(jù)。因此，該階段氣動(dòng)彈性分析以定性為主，通過(guò)建立多學(xué)科綜合分析與優(yōu)化設(shè)計(jì)平臺(tái)，開(kāi)展關(guān)鍵參數(shù)敏感性分析，找出結(jié)構(gòu)參數(shù)優(yōu)化設(shè)計(jì)方向，可用于啟發(fā)或指導(dǎo)后續(xù)詳細(xì)設(shè)計(jì)工作。但在飛機(jī)設(shè)計(jì)方案階段，如何針對(duì)大量不確定參數(shù)開(kāi)展高效分析工作，以加快總體方案迭代進(jìn)度，一直是困擾氣動(dòng)彈性設(shè)計(jì)師的難題[4-5]。

目前工程上氣動(dòng)彈性穩(wěn)定性分析主要針對(duì)顫振問(wèn)題，常用方法為模態(tài)法，即將分析模型轉(zhuǎn)換到模態(tài)空間下求解，因此，模態(tài)分析是顫振分析的必要流程[6-8]。本文提出一種基于虛擬質(zhì)量的飛機(jī)顫振快速分析方法，可減少模態(tài)分析次數(shù)，有效提升分析效率，滿足方案設(shè)計(jì)階段快速分析需求。

模態(tài)法采用一定數(shù)量低頻固有模態(tài)表征結(jié)構(gòu)位移響應(yīng)，由于模態(tài)振型和頻率都處于廣義坐標(biāo)系，與詳細(xì)有限元模型無(wú)關(guān)，因此，可有效降低求解對(duì)象階數(shù)，同時(shí)不會(huì)顯著影響模型主要響應(yīng)特性[9-10]。但是如果求解對(duì)象需要包含局部變形信息，選取的模態(tài)不能表征時(shí)，就可能導(dǎo)致計(jì)算錯(cuò)誤。針對(duì)這類問(wèn)題，虛擬質(zhì)量法提供了一種解決辦法，通過(guò)附加質(zhì)量，獲得局部參數(shù)變化引起的變形信息，提高分析準(zhǔn)確性。虛擬質(zhì)量法最早由Karpel 在子結(jié)構(gòu)模態(tài)綜合中提出[11]，后來(lái)推廣到其他問(wèn)題。如在“濕模態(tài)”求解中，虛擬質(zhì)量可用以體現(xiàn)不可壓縮流體對(duì)結(jié)構(gòu)產(chǎn)生的作用，有效減少工作量[12-13]；在損傷識(shí)別中，基于虛擬質(zhì)量法，有效增加測(cè)試數(shù)據(jù)量，避免結(jié)構(gòu)附加真實(shí)物理元件的困難[14-15]。

1 理論基礎(chǔ)

結(jié)構(gòu)自由度下，氣動(dòng)彈性運(yùn)動(dòng)方程可表示為

虛擬質(zhì)量法是在原結(jié)構(gòu)中施加一塊集中質(zhì)量ΔM[16]，原結(jié)構(gòu)去掉阻尼項(xiàng)的自由振動(dòng)方程變?yōu)?/p>

通過(guò)模態(tài)分析，式（3）可表示為

用虛擬質(zhì)量模態(tài)振型ΦFM替代原結(jié)構(gòu)模態(tài)振型Φ，同時(shí)去掉虛擬質(zhì)量，對(duì)質(zhì)量矩陣、阻尼矩陣、剛度矩陣以及氣動(dòng)力矩陣進(jìn)行廣義化，自由振動(dòng)方程變?yōu)?/p>

顫振方程變?yōu)?/p>

通過(guò)求解式（5）和式（6）可分別得到結(jié)構(gòu)固有模態(tài)頻率和顫振特性。

傳統(tǒng)顫振分析流程需要先進(jìn)行模態(tài)分析，然后，在此基礎(chǔ)上開(kāi)展顫振分析，對(duì)每一組結(jié)構(gòu)參數(shù)都需要重復(fù)這個(gè)過(guò)程。按照上述方程推導(dǎo)，本文方法通過(guò)在原結(jié)構(gòu)基準(zhǔn)模型適當(dāng)位置，施加集中質(zhì)量，開(kāi)展模態(tài)分析得到虛擬質(zhì)量模態(tài)振型，將該振型作為固定的統(tǒng)一振型。虛擬質(zhì)量施加位置一般選取靠近目標(biāo)參數(shù)研究區(qū)域，同時(shí)不產(chǎn)生低頻局部模態(tài)，大小最好與模型總質(zhì)量在一個(gè)數(shù)量級(jí)。當(dāng)結(jié)構(gòu)參數(shù)局部改變時(shí)，不再開(kāi)展模態(tài)分析，直接利用虛擬質(zhì)量模態(tài)振型對(duì)模型質(zhì)量矩陣、阻尼矩陣、剛度矩陣以及氣動(dòng)力矩陣進(jìn)行廣義化，進(jìn)而求解顫振方程。本文方法在參數(shù)研究時(shí)只需開(kāi)展一次模態(tài)分析，節(jié)省了模態(tài)分析時(shí)間，有效提升了顫振分析效率。

2 分析流程

虛擬質(zhì)量法的基本思想是振型統(tǒng)一[17-18]。分析流程的第一步是在基準(zhǔn)有限元模型上施加一個(gè)集中質(zhì)量，作為虛擬質(zhì)量，在此基礎(chǔ)上開(kāi)展模態(tài)分析，得到虛擬質(zhì)量模態(tài)振型；然后，針對(duì)不同參數(shù)下的質(zhì)量矩陣、阻尼矩陣、剛度矩陣以及氣動(dòng)力矩陣，采用虛擬質(zhì)量模態(tài)振型進(jìn)行廣義化；最后，將廣義化后的質(zhì)量矩陣、阻尼矩陣、剛度矩陣以及氣動(dòng)力矩陣代入顫振方程，進(jìn)而求解顫振特性。當(dāng)結(jié)構(gòu)局部參數(shù)發(fā)生小范圍變化時(shí)，不需要再進(jìn)行模態(tài)分析，本文算例變化范圍為50%以內(nèi)，直接采用該虛擬質(zhì)量模態(tài)振型，更新質(zhì)量矩陣、阻尼矩陣、剛度矩陣以及氣動(dòng)力矩陣即可。整個(gè)分析流程中只需要進(jìn)行一次模態(tài)分析，針對(duì)不同結(jié)構(gòu)參數(shù)，可直接開(kāi)展顫振分析，具體流程見(jiàn)圖1。

圖1 分析流程圖Fig.1 Analysis flow

3 算例仿真

3.1 模型基本參數(shù)

以全動(dòng)翼面為對(duì)象，動(dòng)力有限元模型如圖2，大軸提供彎曲剛度，線彈簧和轉(zhuǎn)軸的組合系統(tǒng)模擬旋轉(zhuǎn)剛度。翼面采用單梁式結(jié)構(gòu)，翼肋、主梁及蒙皮均采用玻璃鋼復(fù)材，梁架網(wǎng)格間填充PMI 硬泡沫，基本參數(shù)如表1 所示。

圖2 有限元模型Fig.2 Finite element model

表1 模型參數(shù)Tab.1 Model parameters

3.2 分析結(jié)果

采用虛擬質(zhì)量為10 kg 的集中質(zhì)量單元，施加位置如圖2 所示。通過(guò)模態(tài)分析，可得到虛擬質(zhì)量模態(tài)振型，如圖3 所示。

圖3 虛擬質(zhì)量模態(tài)振型圖Fig.3 Model shape of the fictitious mass

以模型旋轉(zhuǎn)剛度為參數(shù)，研究旋轉(zhuǎn)剛度對(duì)模型顫振特性的影響。旋轉(zhuǎn)剛度由線彈簧和轉(zhuǎn)軸的組合系統(tǒng)模擬，線彈簧基準(zhǔn)值取1 000 N/mm，在此基礎(chǔ)上旋轉(zhuǎn)剛度分別變化20%和50%，剛度參數(shù)研究范圍如表2 所示。

表2 旋轉(zhuǎn)剛度參數(shù)變化范圍Tab.2 Variation range of rotational stiffness

工程中結(jié)構(gòu)模態(tài)分析方法一般采用MSC 公司的NASTRAN 軟件[19-20]，針對(duì)表2 所示5 個(gè)剛度參數(shù)，分別采用NASTRAN 和本文方法得到前5 階結(jié)構(gòu)固有頻率結(jié)果。分析結(jié)果表明，兩種方法固有模態(tài)頻率結(jié)果吻合非常好，最大誤差不超過(guò)4%，見(jiàn)表3。以旋轉(zhuǎn)剛度為參數(shù)開(kāi)展研究，參數(shù)在50%～150%變化，虛擬質(zhì)量模態(tài)振型可有效替代剛度變化對(duì)模態(tài)振型的影響，本文方法計(jì)算得到的結(jié)構(gòu)固有頻率與NASTRAN 計(jì)算結(jié)果基本一致。

表3 固有模態(tài)頻率Tab.3 Natural mode frequency

針對(duì)旋轉(zhuǎn)剛度的參數(shù)變化，進(jìn)一步比較傳統(tǒng)方法和本文方法的顫振特性計(jì)算結(jié)果。亞聲速非定常氣動(dòng)力計(jì)算都使用ZONA6，采用非匹配點(diǎn)顫振分析，固定氣流密度和馬赫數(shù)，使用g 法求解顫振方程[21]，顫振計(jì)算采用前5 階模態(tài)，計(jì)算馬赫數(shù)為0.1，高度為0。不同旋轉(zhuǎn)剛度計(jì)算結(jié)果表明，兩種方法顫振特性計(jì)算結(jié)果基本一致，顫振速度最大誤差為1.7%，顫振頻率最大誤差為2.1%，顫振速度都隨著旋轉(zhuǎn)剛度的增大而升高，見(jiàn)表4。可見(jiàn)本文方法計(jì)算精度可滿足工程分析需求。

表4 顫振特性Tab.4 Flutter characteristics

以彈簧編號(hào)K1對(duì)應(yīng)的顫振分析結(jié)果為例，對(duì)比兩種方法得到的風(fēng)速-阻尼圖和風(fēng)速-頻率圖，如圖4～圖7所示，兩種方法得到的結(jié)構(gòu)阻尼和頻率隨風(fēng)速變化趨勢(shì)基本一致，阻尼曲線由負(fù)到正的穿越都是第3 階模態(tài)，其0 和2%阻尼值對(duì)應(yīng)風(fēng)速值相差非常小，頻率曲線隨風(fēng)速的變化都表現(xiàn)為第1 階模態(tài)和第3 階模態(tài)的靠近，顫振耦合形式都是彎扭耦合。由此可見(jiàn)，本文方法得到的風(fēng)速-阻尼、風(fēng)速-頻率圖與傳統(tǒng)方法基本一致，顫振特性相當(dāng)，本文方法滿足工程研制中顫振分析需求。

圖4 傳統(tǒng)方法風(fēng)速-頻率圖Fig.4 Velocity vs.frequency graph of traditional method

圖5 傳統(tǒng)方法風(fēng)速-阻尼圖Fig.5 Velocity vs.damping graph of traditional method

圖6 本文方法風(fēng)速-頻率圖Fig.6 Velocity vs.frequency graph of method method in this paper

圖7 本文方法風(fēng)速-阻尼圖Fig.7 Velocity vs.damping graph of method method in this paper

4 結(jié)論

（1）結(jié)合工程需求，提出了一種基于虛擬質(zhì)量的飛機(jī)顫振快速分析方法。該方法通過(guò)在原結(jié)構(gòu)基準(zhǔn)模型適當(dāng)位置，施加一定大小的集中質(zhì)量，開(kāi)展模態(tài)分析得到虛擬質(zhì)量模態(tài)振型，采用該振型作為統(tǒng)一振型，針對(duì)不同參數(shù)下的質(zhì)量矩陣、阻尼矩陣、剛度矩陣以及氣動(dòng)力矩陣進(jìn)行廣義化，后續(xù)變參過(guò)程不再開(kāi)展模態(tài)分析，直接求解顫振特性。

（2）在飛機(jī)研制初期參數(shù)敏感性研究過(guò)程中，基于虛擬質(zhì)量的飛機(jī)顫振快速分析方法，采用虛擬質(zhì)量模態(tài)振型，使得顫振方程的模態(tài)振型固定不變，顯著提高了分析效率，并滿足工程分析計(jì)算精度需求，為后續(xù)飛機(jī)詳細(xì)設(shè)計(jì)和改進(jìn)優(yōu)化等工作提供了技術(shù)思路和解決方案。