姜萌 張美艷 郭其威 唐國(guó)安?
(1.復(fù)旦大學(xué) 航空航天系, 上海 200433) (2.上海宇航系統(tǒng)工程研究所, 上海 201108)
轉(zhuǎn)子的振動(dòng)特性長(zhǎng)期受到學(xué)者關(guān)注,從最初經(jīng)典的兩端簡(jiǎn)支、帶剛性圓盤彈性轉(zhuǎn)軸的Jeffcott轉(zhuǎn)子[1],到單盤轉(zhuǎn)子[2],再到含多個(gè)輪盤的多盤轉(zhuǎn)子[3],研究理論和分析方法不斷豐富完善.張文[4]指出,高轉(zhuǎn)速轉(zhuǎn)子中應(yīng)將葉盤看作彈性體,建立柔盤-柔軸的全彈性轉(zhuǎn)子模型進(jìn)行處理.但是,如果葉片-輪盤-軸各部件都采用有限元模型描述,那么完整的轉(zhuǎn)子模型自由度將十分龐大,計(jì)算中須耗費(fèi)大量時(shí)間,在結(jié)構(gòu)優(yōu)化等應(yīng)用中存在困難.因此,對(duì)葉片-輪盤這一轉(zhuǎn)子中的基本構(gòu)件進(jìn)行模型降階依然具有重要意義,通過(guò)降低單級(jí)彈性葉-盤模型的自由度可有效提高整體計(jì)算的效率.
有限單元法是一種普遍應(yīng)用的數(shù)值方法,也成為許多商用計(jì)算程序分析葉盤動(dòng)力學(xué)特性的基本算法.Hsieh和Abel[5]綜合考慮旋轉(zhuǎn)葉-盤的轉(zhuǎn)動(dòng)非線性,分別利用分布和集中質(zhì)量方法建立有限元模型.Genta和Tonoli[6]將葉片組的陀螺效應(yīng)和離心載荷考慮在內(nèi),建立葉盤組合體的有限元模型.周傳月等[7]將Nastran應(yīng)用于某型葉盤組合體的固有振動(dòng)特性分析中.
為縮減模型總自由度、提高計(jì)算效率,一些基于有限元模型的簡(jiǎn)化數(shù)值算法也日益成熟.常見的簡(jiǎn)化算法有三類:利用葉盤回轉(zhuǎn)對(duì)稱特性的矩陣變換[8,9]、假設(shè)模態(tài)法[10]和模態(tài)綜合法[11,12].有學(xué)者則根據(jù)葉-盤結(jié)構(gòu)的特點(diǎn),采用等效的少量自由度模型進(jìn)行建模,提高計(jì)算效率.如Kaza和Kielb[13]等將葉片等效為彈性梁模型,Turhan和Bulut[14]建立了以歐拉伯努利梁來(lái)代替葉片的理想模型.
在優(yōu)化設(shè)計(jì)等應(yīng)用中,含多級(jí)葉盤的轉(zhuǎn)子以及定子、機(jī)匣等組件的發(fā)動(dòng)機(jī)系統(tǒng)模型的規(guī)模龐大,需進(jìn)一步發(fā)展面向多級(jí)葉片-輪盤系統(tǒng)的動(dòng)力學(xué)分析的方法.一些學(xué)者提出了相應(yīng)的分析模型和算法[15,16].Pan等[17]建立了葉盤軸結(jié)構(gòu)混合維度的有限元模型,實(shí)現(xiàn)三維實(shí)體葉片與二維軸對(duì)稱盤軸的耦合模態(tài)分析,能夠較大地提高葉輪類結(jié)構(gòu)的分析效率而不失分析精度.本文作者根據(jù)周向環(huán)繞的單級(jí)葉片組的物理參數(shù)和模態(tài)分析結(jié)果,將葉片組等效為固有特性相近的正交異性環(huán)形板[18],用軸對(duì)稱模型表示單級(jí)葉盤組合體,從而降低模型的自由度.
本文基于之前的研究,提出了新的正交各向異性軸對(duì)稱葉盤組合體模型,能表現(xiàn)出組合體軸向彎曲與周向扭轉(zhuǎn)的振動(dòng)耦合特性.采用的方法是通過(guò)多個(gè)工況條件下的靜力學(xué)等效估算正交異性彈性參數(shù)的初值,再以固有頻率的誤差作為目標(biāo)函數(shù),由極小化過(guò)程確定彈性參數(shù)的終值.
圖1示意的是具體的研究對(duì)象,一個(gè)帶葉冠的葉盤組合體,外圍葉冠和中部輪盤均可以作為回轉(zhuǎn)體建模,不需縮聚.周向排列的多條葉片將被等效成正交異性環(huán)形板,如圖2所示.
圖1 研究對(duì)象:帶冠葉盤組合體模型Fig.1 Object:Bladed disk model with shroud
圖2 葉片組等效前后模型Fig.2 Initial and equivalent model of blade group
(1)
其中,h0和h1是環(huán)板內(nèi)外緣的待定厚度值.
(2)
其中:
(3)
由于初始帶冠葉片在受到軸向力的作用時(shí),不僅會(huì)產(chǎn)生軸向位移,同時(shí)發(fā)生周向扭轉(zhuǎn).在估算環(huán)板彈性參數(shù)前,需對(duì)其選定一個(gè)材料主軸坐標(biāo)系,使得等效環(huán)板在受相同外力作用時(shí),位移的大小和方向能夠保持和初始葉片一致.為簡(jiǎn)便起見,文中采用單一的角度,在總體上模擬出軸向的加載可產(chǎn)生周向位移的情形.
如圖3所示,分別將各向異性環(huán)板的材料主軸坐標(biāo)系和總體坐標(biāo)系記為o1-x1y1z1、o-xyz,其各個(gè)子午面上材料的主軸y1總是指向環(huán)板的徑向,而x1與環(huán)板所在平面法線x有一個(gè)夾角α.可對(duì)初始帶冠葉片和等效環(huán)板施加相同的外力,通過(guò)調(diào)節(jié)α的大小,保證二者位移效果的一致性,以此確定α和材料主軸坐標(biāo)系.
圖3 環(huán)板材料主軸方向示意圖Fig.3 Material principle axis diagram of annular plate
在材料主軸坐標(biāo)系下,等效環(huán)板的正交各向異性本構(gòu)關(guān)系為:
σ=D0ε
(4)
其中,彈性矩陣D0為:
(5)
含有d11,d12,…,d66共9個(gè)待定的彈性參數(shù).利用轉(zhuǎn)軸公式[19]可將上述主軸坐標(biāo)系下的彈性矩陣轉(zhuǎn)換到總體坐標(biāo)系下:
D=TσD0TσT
(6)
其中,Tσ為應(yīng)力轉(zhuǎn)換矩陣,與夾角α有關(guān).
s.t.x∈X
(7)
其中,X為彈性參數(shù)的可行域集合,C為振動(dòng)模態(tài)的階次集合.
圖4 等效帶冠彈性環(huán)板模型Fig.4 Equivalent elastic annular plate model with shroud
對(duì)于非線性規(guī)劃問(wèn)題(7),目標(biāo)函數(shù)f(x)是非線性的,可能存在多個(gè)局部解,通常需要給出合理的設(shè)計(jì)變量初值,才能獲得有效的極小解.第2節(jié)將給出彈性參數(shù)初值的估算方法.
模型等效過(guò)程示意如圖5所示,其中,深灰色實(shí)體為等效前的周期重復(fù)排列的葉片組,淺灰色實(shí)線圍成的區(qū)域(覆蓋葉片組及其間隙的六面體)為等效后的連續(xù)彈性板.材料主軸o1-x1y1z1由總體坐標(biāo)系o-xyz繞y軸旋轉(zhuǎn)α得到.
圖5 模型等效過(guò)程示意圖Fig.5 Diagram of equivalent model
確定彈性參數(shù)初值時(shí),不必太精確,僅為給出式(7)中設(shè)計(jì)變量的合理區(qū)間.根據(jù)初始葉片組的實(shí)際分布和葉片間的傳力特征,如葉片沿方向o1z1獨(dú)立分布,o1z1方向的正應(yīng)變不能產(chǎn)生該方向的正應(yīng)力,可假設(shè)式(5)彈性矩陣D0中的d33=0,同理,假設(shè)d55=d23=d13=0.考慮到d11與三維板的橫向振動(dòng)有關(guān),在后續(xù)0、1節(jié)徑的模態(tài)中影響不顯著,優(yōu)化時(shí),可暫將初值假設(shè)為:
d11≈d22
(8)
d12是x方向應(yīng)變與其引起的y方向的應(yīng)力之比,改變d12的大小,對(duì)最終頻率的影響亦不顯著,優(yōu)化前可將其初值設(shè)為:
(9)
故僅剩三個(gè)待定的彈性參數(shù)d22、d44和d66,可通過(guò)不同工況下的靜力學(xué)等效估算其初值,靜力學(xué)等效原則為:等效前后的模型受相同外力,產(chǎn)生相同的變形.算例中通過(guò)施加外力fx和外力矩mz,依據(jù)上述等效原則可求得初值,等效時(shí)基于單條葉片,實(shí)現(xiàn)過(guò)程基于最小勢(shì)能原理.
公式推導(dǎo)中,可將回轉(zhuǎn)排列的葉片組近似看作平行排列,那么單條葉片的等效及彈性參數(shù)初值的計(jì)算公式推導(dǎo)均可在直角坐標(biāo)系下完成.將等效單條葉片的彈性板沿周向復(fù)制若干次即可得到整體的等效環(huán)板模型.
假定等效葉片組的三維板,在受到一定外力(包含彎矩、拉力、剪力)作用下的變形位移u是關(guān)于未知系數(shù)a的線性函數(shù),關(guān)系式為:
u=N(x,y,z)a
(10)
其中:
u={u,v,w}T,a={a1,a2,a3}T
N(x,y,z)的形式可根據(jù)初始葉片受力后的實(shí)際變形設(shè)定.根據(jù)幾何方程,得應(yīng)變場(chǎng):
ε=Pa
(11)
利用本構(gòu)方程可得應(yīng)力場(chǎng),由于剛度矩陣中包含待定的彈性參數(shù),故將應(yīng)力場(chǎng)表達(dá)為關(guān)于彈性參數(shù)的函數(shù)形式,如下:
σ(d11,d12,…,d66)=Dε
(12)
積分可得到彈性應(yīng)變能:
(13)
體積分域Ω根據(jù)等效三維板的實(shí)際體積決定.將式(11)、(12)代入式(13)中,可進(jìn)一步將彈性應(yīng)變能整理成關(guān)于未知系數(shù)a的二次型形式:
(14)
其中:
(15)
若在與葉片與葉冠相連的面上施加外力f,產(chǎn)生位移u,則外力對(duì)等效三維板的做功大小為:
(16)
其中:
(17)
面積分域S由等效三維板受力面的實(shí)際大小決定.若在與葉片與葉冠相連的面上施加外力矩m,產(chǎn)生轉(zhuǎn)角θ,則在計(jì)算做功時(shí),分別用m,θ替換式(16)中的f,u即可.
綜上,可得環(huán)板受外力作用變形后的總勢(shì)能表達(dá)式如下:
(18)
依據(jù)最小勢(shì)能原理的駐值條件[19],可求得:
a(d11,d12,…,d66)=A-1b
(19)
將式(19)代入(10)中,得到僅含有未知彈性參數(shù)的位移表達(dá)式:
u(d11,d12,…,d66)=N(x,y,z)A-1b
(20)
(21)
即可構(gòu)造方程組,求解后得到待定彈性參數(shù)值,以此作為三維環(huán)板彈性參數(shù)的初值.
算例的初始葉片組模型如圖1,共含有72條葉片,每隔5°均勻分布.等效時(shí),將每隔5°分布的單條葉片等效為占據(jù)5°的彈性板,如圖6中所示環(huán)板上標(biāo)注的5°扇區(qū)陰影區(qū)域.完成單個(gè)扇區(qū)的等效建模后,繞軸旋轉(zhuǎn)復(fù)制71次即得到葉片組等效后的三維環(huán)板模型.
圖6 等效葉片組的三維帶冠環(huán)板幾何尺寸示意圖Fig.6 Geometric size diagram of three-dimensional annular plate of equivalent blade group
初始葉片組模型的幾何與材料等物理參數(shù)值見表1,均采用國(guó)際單位制.
表1 葉片組模型物理參數(shù)表Table 1 Parameters of blade group
葉冠為各向同性的軸對(duì)稱柱殼,幾何與材料等物理參數(shù)值見表2,模型不需等效.
表2 葉冠模型參數(shù)表Table 2 Parameters of shroud
等效單條葉片的三維彈性板,密度ρ取為與初始葉片組相等,即ρ=ρb,y方向(徑向)的長(zhǎng)度Y為葉片內(nèi)、外徑之差,外緣z方向的長(zhǎng)度Z可用弧長(zhǎng)近似,即:
Y=R1-R0=0.1191m
Z≈R1θ=0.031m
(22)
式中θ=π/36,將表1數(shù)據(jù)代入公式(1)~(3)中,經(jīng)計(jì)算可得等效環(huán)板的厚度分布函數(shù)為:
(23)
基于第2節(jié)分析思路,首先確定材料坐標(biāo)系與總體坐標(biāo)系間的夾角α=66°.將α=66°代入轉(zhuǎn)軸公式中,計(jì)算得總體坐標(biāo)系下的剛度矩陣D為:
(24)
其中:
以下分別施加外力矩mz和外力fx,通過(guò)靜力學(xué)等效求得d22、d44和d66的初值.
第一步,構(gòu)造純彎曲變形,確定d22.葉根固定,葉冠處受沿z方向的力矩:
(25)
靜力分析后,初始單條葉片在xoy和yoz平面的變形如圖7所示.黑色虛線為變形前葉片組位置,藍(lán)色實(shí)線為變形后葉片組位置.可見,葉片x方向幾乎無(wú)變形位移,在yoz平面的變形可近似看成歐拉梁,其中z方向的變形可用二次函數(shù)模擬,因此可設(shè)等效三維板的位移場(chǎng)為:
(26)
將式(25)、(26)代入(10)~(20)可得等效板的位移場(chǎng)表達(dá)式,其中包含了待定彈性參數(shù)d22.
(27)
將等效環(huán)板最大y坐標(biāo)所在平面上的中心點(diǎn)坐標(biāo)和位移數(shù)值代入(27),可求解出d22,結(jié)果為:
d22=2.9×1012Pa
(28)
圖7 受外力矩作用的葉片xoy和yoz平面變形圖Fig.7 Deformation diagram of blade under moment in the plane of xoy and yoz
第二步,構(gòu)造純剪變形,確定d44和d66.葉根固定,葉冠處受沿x方向的外力:
(29)
靜力分析后,初始單條葉片在xoy和yoz平面的變形如圖8所示.黑色虛線為變形前葉片組位置,藍(lán)色實(shí)線為變形后葉片組位置.
可見,葉片在y方向無(wú)位移,x,z方向的變形可用三次函數(shù)模擬,故設(shè)等效板的位移場(chǎng)為:
(30)
將式(29)、(30)代入式(10)~(20)可得關(guān)于d44,d66的位移表達(dá)式,將等效環(huán)板最大y坐標(biāo)所在平面上的中心點(diǎn)坐標(biāo)和位移數(shù)值代入可求得:
d44=6.5×1014Pa
d66=1.9×107Pa
(31)
圖8 受外力矩作用的葉片xoy和yoz平面變形圖Fig.8 Deformation diagram of blade under moment in the plane of xoy and yoz
參見第2.1節(jié)式(8)、(9),可將d11,d12的初值取為:
d11=d22=2.9×1012Pa
d12=1.1×1012Pa
(32)
根據(jù)式(7)、設(shè)置設(shè)計(jì)變量的可行域,建立數(shù)學(xué)規(guī)劃模型.
s.t.d∈[106,1015]×[106,1015]×
[106,1015]×[106,1015]×[106,1015]
(33)
其中,d={d11,d12,d22,d44,d66},將式(28)、(31)和(32)中的數(shù)據(jù)作為初值,利用MSC.Nastran的動(dòng)力優(yōu)化功能[20]得彈性參數(shù)終值.迭代計(jì)算后可確定環(huán)板彈性參數(shù)的終值為:
d11=d22=4.1×1012Pa
d12=3.5×1010Pa
d44=8.2×1014Pa
d66=6.1×107Pa
(34)
利用上述方法求得三維環(huán)板建模所需的全部幾何、材料參數(shù)后,可在Patran中建立環(huán)板的有限元模型,將等效三維環(huán)板模型與初始輪盤組合后可得到三維葉盤縮聚模型.
模態(tài)分析后得到的葉盤縮聚模型與初始葉盤模型的0、1節(jié)徑模態(tài)變形云圖如圖9,可見模態(tài)變形的橫向等位移線基本一致.
圖9 初始與三維葉盤縮聚模型模態(tài)云圖對(duì)比Fig.9 Comparison of modal deformation nephogram for initial and reduced bladed disk
利用兩種模型計(jì)算得到的固有頻率結(jié)果如表3,其中0節(jié)徑和1節(jié)徑的固有頻率相差分別為-4.18%和1.31%.算例中,初始葉盤模型為回轉(zhuǎn)周期結(jié)構(gòu),模態(tài)分析需用到單個(gè)扇區(qū)的模型,節(jié)點(diǎn)數(shù)目約為3100;等效后的縮聚模型為軸對(duì)稱分布,模態(tài)分析時(shí)僅需用到子午面的模型,節(jié)點(diǎn)數(shù)目不足350,可大大降低模型動(dòng)力學(xué)分析時(shí)的自由度數(shù)目,達(dá)到高效計(jì)算的目的.
表3 初始與三維板縮聚模型固有頻率數(shù)據(jù)結(jié)果Table 3 Natural frequencies of initial and reduced model
文章[18]建立的二維葉盤縮聚模型與本文建立三維葉盤縮聚模型模態(tài)變形后的單元輪廓對(duì)比圖如圖10所示,黑色線條和藍(lán)色線條分別為模態(tài)變形前、后的單元輪廓線,可見,三維葉盤縮聚模型與二維葉盤縮聚模型相比,更能正確表現(xiàn)葉片-輪盤的軸向-周向振動(dòng)耦合效應(yīng).分析表明,文中的等效建模方法能相當(dāng)可觀的降低葉盤整體有限元分析的自由度數(shù)目,且用縮聚模型計(jì)算得到的低階振動(dòng)固有頻率具有的足夠精度.
圖10 二維與三維葉盤縮聚模型模態(tài)變形單元輪廓圖Fig.10 Comparisonof modal deformation nephogram for 2D & 3D reduced bladed disk
在較小的尺度上,葉片組是由葉片按照固定回轉(zhuǎn)角度沿周向?qū)ΨQ排列的周期性重復(fù)結(jié)構(gòu).在較大的尺度上,這種周期性重復(fù)結(jié)構(gòu)可以被近似地均質(zhì)化,這樣可使得分析和計(jì)算難度大大降低.對(duì)于葉片組的等效建模思路,本文選用各項(xiàng)異性環(huán)形板等效周期性重復(fù)的葉片組三維有限元模型,在保證兩種模型外觀特征尺寸一致、整體質(zhì)量和慣量相同的條件下,通過(guò)調(diào)節(jié)等效環(huán)板模型的彈性模量等參數(shù),使得初始模型和等效模型的主要固有頻率和模態(tài)具有相似的特征.
實(shí)現(xiàn)模型的等效時(shí),綜合運(yùn)用了彈性力學(xué)中的最小勢(shì)能等原理,給出了具有可操作性的確定等效模型參數(shù)的方法,并以航空發(fā)動(dòng)機(jī)帶葉冠的低壓渦輪葉片-輪盤模型為具體應(yīng)用實(shí)例,展示了模型等效的步驟.對(duì)于設(shè)計(jì)狀態(tài)已確定的單級(jí)葉盤,等效的環(huán)板和輪盤組合體具有軸對(duì)稱特點(diǎn),故計(jì)算效率高,可直接應(yīng)用于轉(zhuǎn)子的整體動(dòng)力學(xué)計(jì)算,也可為轉(zhuǎn)子系統(tǒng)的優(yōu)化設(shè)計(jì)等應(yīng)用提供精度和效率兼顧的彈性葉盤縮聚模型.