石先杰， 左 朋,2

(1. 中國工程物理研究院 總體工程研究所，四川 綿陽 621999；2. 中國科學(xué)技術(shù)大學(xué) 近代力學(xué)系，合肥 230026)

功能梯度材料[1-3]作為一種新興材料，它克服了傳統(tǒng)材料在過熱環(huán)境下機(jī)械性能下降的缺陷，具有極強(qiáng)的耐熱性，而由其所制成的環(huán)板結(jié)構(gòu)有著質(zhì)量輕、結(jié)構(gòu)剛度大、力學(xué)性能好、承載能力強(qiáng)等優(yōu)點(diǎn)，在航空航天工程、飛機(jī)工程、核反應(yīng)堆等領(lǐng)域的發(fā)動(dòng)機(jī)端蓋上有著廣泛應(yīng)用。因此，這類結(jié)構(gòu)在實(shí)際工程應(yīng)用中通常會(huì)受到復(fù)雜工作環(huán)境的影響而產(chǎn)生較大的振動(dòng)甚至引起疲勞損傷。因此，準(zhǔn)確分析預(yù)測(cè)功能梯度環(huán)板結(jié)構(gòu)振動(dòng)特性一直是國內(nèi)外專家學(xué)者的研究熱點(diǎn)。為此，國內(nèi)外學(xué)者已開展了一系列數(shù)值求解方法研究，如Rayleigh-Ritz法[4-6]、微分求積法[7-10]、擬格林函數(shù)法[11]和有限環(huán)形棱柱法[12]等。

然而，上述研究均是基于常溫環(huán)境，在高溫環(huán)境下，結(jié)構(gòu)的材料屬性會(huì)隨著溫度的變化而改變，從而導(dǎo)致其自由振動(dòng)的固有頻率相比于常溫下變得明顯不同。對(duì)于熱環(huán)境下功能梯度結(jié)構(gòu)的振動(dòng)特性，國內(nèi)外學(xué)者進(jìn)行了一些研究分析。假設(shè)材料特性依賴于溫度，并在板的厚度方向上逐漸變化，Haddadpour等[13]采用伽遼金法求解運(yùn)動(dòng)方程，對(duì)簡支功能梯度圓柱殼在熱環(huán)境中的自由振動(dòng)特性進(jìn)行了分析。Mirtalaie[14]根據(jù)微分求積法，求解了熱環(huán)境下功能梯度薄環(huán)形扇形板自由振動(dòng)頻率。在三維彈性理論的框架下，Yang等[15]分析了功能梯度環(huán)板的熱振響應(yīng)。在一階剪切變形理論的框架內(nèi)應(yīng)用微分求積法，Malekzadeh等[16-17]研究了旋轉(zhuǎn)功能梯度圓柱殼和帶有彈性支撐的功能梯度環(huán)板在熱場(chǎng)中的自由振動(dòng)特性。用傅里葉級(jí)數(shù)結(jié)合輔助函數(shù)來描述板的位移容許函數(shù)，Zhou等[18]研究了彈性地基上具有溫度梯度的超音速多孔功能梯度材料板的振動(dòng)和顫振特性。Li等[19]采用特征正交多項(xiàng)式詳細(xì)研究了功能梯度多孔階梯圓柱殼在熱場(chǎng)中的自由振動(dòng)、穩(wěn)態(tài)和瞬態(tài)響應(yīng)特性。蒲育等[20-21]根據(jù)二維彈性理論，采用微分求積法分析了材料屬性徑向變化的功能梯度薄環(huán)板在熱環(huán)境下的面內(nèi)自由振動(dòng)特性。同樣是考慮材料屬性沿徑向變化功能梯度薄環(huán)板，呂朋等[22]采用改進(jìn)傅里葉級(jí)數(shù)建立了其在熱場(chǎng)中受彈性約束的面內(nèi)振動(dòng)分析模型。而對(duì)于沿厚度方向變化的功能梯度環(huán)板，石鵬等[23]應(yīng)用三維Chebyshev-Ritz法對(duì)其在熱環(huán)境下的自由振動(dòng)特性進(jìn)行了求解。

綜上所述，目前關(guān)于功能梯度環(huán)板在熱場(chǎng)中的自由振動(dòng)分析雖然有一定的研究成果，但大部分研究工作集中在材料屬性沿徑向變化的薄環(huán)板面內(nèi)振動(dòng)問題。而對(duì)于厚環(huán)板問題，往往是采用三維彈性理論求解，存在著計(jì)算效率低等問題。因此，針對(duì)工程中更為常見的中厚環(huán)板，考慮材料特性沿厚度方向的梯度變化，本文以一階剪切變形理論為基礎(chǔ)，通過沿邊界設(shè)置具有可變剛度的邊界約束彈簧的方式來考慮任意邊界約束。同時(shí)，采用譜幾何法建立功能梯度環(huán)板熱振分析模型。在譜幾何法求解框架下，板的位移函數(shù)總是表示為標(biāo)準(zhǔn)的傅里葉余弦級(jí)數(shù)，并輔以四個(gè)正弦函數(shù)項(xiàng)來消除邊界處位移函數(shù)的導(dǎo)數(shù)不連續(xù)性問題。這種級(jí)數(shù)展開方法理論上能夠表征任意有限長結(jié)構(gòu)的位移變量，且函數(shù)形式統(tǒng)一，易于參數(shù)化分析。在數(shù)值算例中，通過與文獻(xiàn)解、有限元解的對(duì)比研究，驗(yàn)證了該模型的正確性，進(jìn)而分析了邊界條件、環(huán)境溫度、梯度指數(shù)等對(duì)功能梯度環(huán)板固有頻率的影響。

1 理論方法

1.1 模型描述

考慮如圖1所示的功能梯度環(huán)板，圓柱坐標(biāo)系(r,θ,z)位于結(jié)構(gòu)的幾何中面上。圖1中：R0和R1分別為環(huán)板沿徑向方向r的內(nèi)徑和外徑；θ為圓周方向；h為環(huán)板的厚度。文中研究的功能梯度材料由陶瓷和金屬材料兩種成分組成，其有效材料參數(shù)P(包括彈性模量E、密度ρ、泊松比ν和熱膨脹系數(shù)α等)沿環(huán)板厚度方向呈梯度變化，可以被描述為[24]

(1)

式中：Pt和Pu分別為環(huán)板上、下表面的材料參數(shù)；-h/2≤z≤h/2；p為功能梯度指數(shù)?？紤]熱環(huán)境的影響，功能梯度材料的有效材料參數(shù)Pξ(ξ=t,u)在溫度值T的變化影響下有著如式(2)所示的性質(zhì)

Pξ(z,T)=P0ξ(P-1ξT-1+1+P1ξT+P2ξT2+P3ξT3)

(2)

式中，P-1ξ,P0ξ,P1ξ,P2ξ和P3ξ為溫度相關(guān)系數(shù)。文中研究的是均勻溫度分布的情形[25]，T0為參考溫度，其溫度變化為ΔT=T-T0。

圖1 功能梯度環(huán)板的幾何形狀以及坐標(biāo)系Fig.1 Geometric shape and coordinate system of functionally graded annular plate

1.2 能量方程

根據(jù)一階剪切變形理論，環(huán)板上任一點(diǎn)沿r,θ和z方向的位移場(chǎng)分量u,v,w的表達(dá)式可以寫成

u(r,θ,z,t)=u0(r,θ,t)+zφr(r,θ,t),
v(r,θ,z,t)=v0(r,θ,t)+zφθ(r,θ,t),
w(r,θ,z,t)=w0(r,θ,t)

(3)

式中：u0,v0,w0分別為環(huán)板在幾何中面上沿r,θ和z方向的位移分量；φr和φθ分別為在rz和θz平面上的轉(zhuǎn)角位移分量。

根據(jù)中厚板理論，環(huán)板結(jié)構(gòu)中的應(yīng)變和位移之間存在線性關(guān)系式

(4)

(5)

ε=[εr,εθ,γrθ,γrz,γθz]T,
σ=[σr,σθ,τrθ,τrz,τθz]T,
σ=Q(5×5)ε

(6)

式中，Q(5×5)為與溫度值T和厚度坐標(biāo)z有關(guān)的彈性常數(shù)所組成的矩陣(參考2011年Malekzadeh等的研究)，其元素的具體值為

(7)

對(duì)應(yīng)力分量σ在z方向積分，引入剪切修正系數(shù)κ=5/6, 得到環(huán)板的本構(gòu)方程表達(dá)式

N=Aε0+Bχ,M=Bε0+Dχ,Qκ=Aκγ

(8)

式中：N=[Nr,Nθ,Nrθ]T和Q=[Qr,Qθ]T分別為環(huán)板的合力和橫向剪應(yīng)力分量；M=[Mr,Mθ,Mrθ]T為合彎矩分量；A(Aκ),B和D分別為以結(jié)構(gòu)的拉伸剛度、拉伸-彎曲耦合剛度和彎曲剛度為元素組成的矩陣，其具體表達(dá)式可見2011年Malekzadeh等的研究。

根據(jù)建立的本構(gòu)方程表達(dá)式，結(jié)合2011年Malekzadeh等給出的熱應(yīng)力描述，得到環(huán)板的應(yīng)變能表達(dá)式

(9)

代入位移場(chǎng)分量可以得到環(huán)板的動(dòng)能表達(dá)式

(10)

文中通過在環(huán)板的內(nèi)部和外部放置邊界彈簧實(shí)現(xiàn)邊界條件的模擬。邊界約束彈簧所儲(chǔ)存的能量可以寫成

(11)

式中，上標(biāo)0和1分別為設(shè)置在環(huán)板內(nèi)部和外部的邊界彈簧剛度值。

1.3 位移容許函數(shù)以及方程的求解

文中分別采用了譜幾何法和傅里葉正余弦函數(shù)來表示環(huán)板沿徑向方向和周向方向的位移容許函數(shù)

(12)

如上所述，熱環(huán)境下功能梯度環(huán)板的拉格朗日能量泛函可以寫成

(13)

根據(jù)1.2節(jié)的能量方程表達(dá)式，結(jié)合相應(yīng)的位移容許函數(shù)進(jìn)行求解，并采用Rayleigh-Ritz法對(duì)未知展開系數(shù)求偏導(dǎo)，可以得到熱環(huán)境下功能梯度環(huán)板自由振動(dòng)特征方程

(K-ω2M)H=0

(14)

式中：K和M分別為環(huán)板的剛度矩陣和質(zhì)量矩陣；H為未知系數(shù)向量，其中溫度變化會(huì)影響K。

2 數(shù)值算例與分析

根據(jù)上述分析獲得的自由振動(dòng)特征方程，本章將分析熱環(huán)境下功能梯度環(huán)板自由振動(dòng)特性。根據(jù)航空航天飛行器的工程實(shí)際應(yīng)用，在接下來的數(shù)值算例中選擇了幾何尺寸為R0=0.5 m,R1=2.0 m,h=0.10 m的環(huán)板結(jié)構(gòu)作為研究對(duì)象。通過設(shè)置邊界約束彈簧的剛度值，文中研究考慮了包括經(jīng)典和彈性邊界在內(nèi)的10種任意邊界條件，并且采用了不同字母來表示，其中： C，SS，SD和F分別為固支、簡支、剪切-膜片邊界和自由邊界； E1，E2，E3，E4，E5和E6為彈性邊界。它們的彈簧剛度值按Su等的研究設(shè)置。此外，位移容許函數(shù)在徑向和周向方向的級(jí)數(shù)在理論上可以展開為無窮項(xiàng)。但是考慮到實(shí)際數(shù)值計(jì)算需要，文中對(duì)其截?cái)嘀劣邢揄?xiàng)。參考文獻(xiàn)[26]對(duì)譜幾何法收斂性分析結(jié)果，文中研究將徑向和周向的級(jí)數(shù)均截?cái)酁镸=N=14。

2.1 數(shù)值驗(yàn)證

不同邊界條件和梯度指數(shù)p下功能梯度環(huán)板在不考慮熱環(huán)境影響情況下的第1階固有頻率，如表1所示。算例包括了兩種板的厚度h=0.05 m和h=0.10 m, 其余幾何尺寸與默認(rèn)參數(shù)保持一致。環(huán)板的上、下表面材料分別為鋁和氧化鋯，材料參數(shù)為:Et=168 GPa,νt=0.3,ρt=5 700 kg/m3,Eu=70 GPa,νu=0.3,ρu=2 707 kg/m3。表1中將譜幾何法模型計(jì)算結(jié)果與Su等的研究結(jié)果進(jìn)行了對(duì)比。從表1可以看出，文中結(jié)果與Su等的解非常接近，存在的偏差主要來源于徑向所采用的位移函數(shù)差異。這表明文中模型邊界約束彈簧剛度以及截?cái)鄶?shù)的設(shè)置是合理的，可以有效求解具有不同梯度指數(shù)p的功能梯度環(huán)板的自由振動(dòng)特性。

在開展熱環(huán)境下功能梯度環(huán)板結(jié)構(gòu)自由振動(dòng)分析時(shí)，由于缺乏相關(guān)文獻(xiàn)數(shù)據(jù)作為參考。因此，文中采用有限元法(finite element method，F(xiàn)EM)計(jì)算結(jié)果作為對(duì)比數(shù)據(jù)。環(huán)板上、下表面的材料分別為氮化硅(Si3N4)和304不銹鋼(SUS304)，它們?cè)跓釄?chǎng)中的材料參數(shù)如表2所示。環(huán)板的幾何尺寸與默認(rèn)參數(shù)保持一致，參考溫度取為T0=300 K。同樣是考慮了C-C，SS-SS和C-SS三種經(jīng)典邊界條件。不同溫度變化ΔT下環(huán)板前5階固有頻率對(duì)比情況(梯度指數(shù)p=1.0)，如表3所示。不同梯度指數(shù)p下的前9階固有頻率對(duì)比(溫度變化ΔT=150 K)，如表4所示。兩種方法得到的頻率結(jié)果可以較好地匹配，最大偏差在0.8%以內(nèi)。偏差的來源主要是文中方法和FEM分別采用了不同的板殼理論和求解方法體系。此外，F(xiàn)EM計(jì)算結(jié)果的求解精度依賴于網(wǎng)格密度，需要不斷地細(xì)化網(wǎng)格才能獲得更加準(zhǔn)確的結(jié)果，但同時(shí)也會(huì)帶來計(jì)算成本的增加。兩種方法求解獲得的模態(tài)振型對(duì)比情況，如圖2所示。其中溫度變化ΔT=150 K, 梯度指數(shù)p=0.5, 邊界條件為C-C。從圖2可以看出，基于兩種方法獲得的模態(tài)振型一致性較好。綜上所述，文中方法可以有效地分析熱環(huán)境下功能梯度環(huán)板自由振動(dòng)特性。

表1 具有不同邊界條件和梯度指數(shù)的功能梯度環(huán)板的第1階固有頻率Tab.1 The first natural frequency of functionally graded annular plate with different boundary conditions and gradient indexes

2.2 參數(shù)化分析

在驗(yàn)證文中構(gòu)建模型準(zhǔn)確性和通用性的基礎(chǔ)上，針對(duì)航空航天飛行器工程運(yùn)用中的各種復(fù)雜工況，本節(jié)對(duì)功能梯度環(huán)板在熱環(huán)境下的振動(dòng)特性開展參數(shù)化分析。結(jié)構(gòu)的幾何尺寸以及材料的選擇與表3、表4保持一致。

首先，當(dāng)溫度變化ΔT=50 K時(shí)，各種復(fù)雜邊界條件下功能梯度環(huán)板固有頻率，如表5所示，材料梯度指數(shù)p=1。從表5可以看出，對(duì)于經(jīng)典邊界，C-SD邊界下的固有頻率最高，SS-SD和SD-SD邊界下的頻率結(jié)果非常接近，而F邊界的存在會(huì)很明顯地降低結(jié)構(gòu)的固有頻率。而對(duì)于彈性邊界，在E1-E1和E4-E4邊界條下環(huán)板分別有著最高和最低的固有頻率。

溫度變化ΔT對(duì)功能梯度環(huán)板固有頻率(n=1,m=1)的影響情況，如圖3所示。圖3中，材料的梯度指數(shù)p=1。不難看出，在經(jīng)典、經(jīng)典-彈性以及彈性邊界條件下，隨著ΔT增大，環(huán)板的固有頻率會(huì)不斷降低。這是因?yàn)棣變大降低應(yīng)變能導(dǎo)致了功能梯度環(huán)板剛度降低，從而使結(jié)構(gòu)固有頻率減小。

梯度指數(shù)p以及周向波數(shù)n的變化對(duì)環(huán)板固有頻率的影響，如圖4所示，溫度變化ΔT=50 K。從圖4可以發(fā)現(xiàn)，結(jié)構(gòu)固有頻率會(huì)隨著周向波數(shù)n的增加而不斷增加。而隨著梯度指數(shù)p的增加，功能梯度環(huán)板的固有頻率先是迅速下降，達(dá)到一定值后，影響逐漸不明顯。因?yàn)殡S著p的增加，環(huán)板材料屬性更加接近于下表面的金屬組分，從而降低了結(jié)構(gòu)的固有頻率。

圖4 周向波數(shù)n和梯度指數(shù)的變化對(duì)功能梯度環(huán)板固有頻率的影響Fig.4 The influence of changes in circumferential wavenumber n and gradient indexes on the natural frequencies of functionally graded annular plates

3 結(jié) 論

綜合考慮功能梯度環(huán)板材料屬性溫度相關(guān)性及其在厚度方向上梯度連續(xù)變化，本文在一階剪切變形理論框架下推導(dǎo)出功能梯度環(huán)板振動(dòng)能量方程，并采用譜幾何法求解得到功能梯度環(huán)板在熱場(chǎng)中的自由振動(dòng)特性。通過將文中計(jì)算結(jié)果與相關(guān)文獻(xiàn)和有限元法結(jié)果進(jìn)行對(duì)比，驗(yàn)證了所構(gòu)建分析模型的正確性和可靠性。在此基礎(chǔ)上，研究了邊界約束條件、溫度變化和梯度指數(shù)等參數(shù)對(duì)熱環(huán)境下功能梯度環(huán)板自由振動(dòng)特性的影響，獲得以下結(jié)論：

(1) 文中所構(gòu)建的計(jì)算分析模型可以有效分析預(yù)示熱環(huán)境下功能梯度環(huán)板自由振動(dòng)特性。

(2) 對(duì)于經(jīng)典邊界，SS和SD邊界對(duì)功能梯度環(huán)板自由振動(dòng)固有頻率的影響類似，F(xiàn)邊界則會(huì)很明顯地降低結(jié)構(gòu)固有頻率。對(duì)于彈性邊界，E1和E4邊界下的功能梯度環(huán)板分別有著最高和最低的固有頻率。

(3) 功能梯度環(huán)板固有頻率隨著溫度變化ΔT的增加而降低；而梯度指數(shù)p的增加，則會(huì)導(dǎo)致功能梯度環(huán)板的固有頻率呈現(xiàn)先迅速下降后逐漸穩(wěn)定的趨勢(shì)。