粘貼壓電層功能梯度材料Euler梁熱屈曲前自由振動(dòng)的半解析數(shù)值法

蘇厚德， 樊建領(lǐng)， 俞樹(shù)榮， 馮玉潔， 曹維國(guó)

(1. 甘肅藍(lán)科石化高新裝備股份有限公司，蘭州730070;2.蘭州理工大學(xué) 石油化工學(xué)院，蘭州730050)

用壓電材料制成的作動(dòng)器在結(jié)構(gòu)智能控制中已經(jīng)具有廣泛的應(yīng)用。關(guān)于壓電智能結(jié)構(gòu)控制的研究工作主要涉及梁、板、殼結(jié)構(gòu)的靜動(dòng)態(tài)響應(yīng)控制［1～17］。Crawley 和Luis［1］對(duì)于壓電層合梁做出有意義工作，給出了壓電材料與結(jié)構(gòu)體之間的靜動(dòng)態(tài)藕合模型，并通過(guò)實(shí)驗(yàn)驗(yàn)證了理論結(jié)果的可靠性，但是在Crawley 的梁模型中沒(méi)有考慮橫向剪切變形以及軸向力的影響，Zhou 等［2］，Xia 等［3］對(duì)表面粘貼壓電片的梁式板和變厚度板結(jié)構(gòu)的幾何非線性大變形及振動(dòng)響應(yīng)智能控制進(jìn)行了研究。丁麗霞等［4］考慮上下嵌有功能梯度原件的薄板，基于經(jīng)典板理論，導(dǎo)出了基于功能梯度彈性薄板小撓度屈曲平衡微分方程，研究了材料的梯度指數(shù)對(duì)于臨界電壓的影響。Dai 等［5］基于一階剪切理論同時(shí)考慮熱-電-機(jī)耦合效應(yīng)，用無(wú)網(wǎng)格Galerkin 法對(duì)含有壓電驅(qū)動(dòng)和傳感的FGM 板靜動(dòng)態(tài)響應(yīng)做了研究。林啟榮等［6］分析了壓電驅(qū)動(dòng)器對(duì)彈性梁的等效作用，并給出了一端固支一端簡(jiǎn)支梁的解析解。Fridman 等［7］基于一階剪切理論，研究了對(duì)稱(chēng)和非對(duì)稱(chēng)鋪設(shè)壓電層的柔性層合復(fù)合梁的固有頻率，屈曲載荷等與軸向、厚度方向載荷及電場(chǎng)載荷之間的特征曲線關(guān)系。Shen等［8，9］、Huang 等［10］基于Reddy 高階剪切理論和von Karman 幾何非線性方程，同時(shí)考慮物性參數(shù)的溫度依賴(lài)效應(yīng)，采用攝動(dòng)法對(duì)壓電-功能梯度板結(jié)構(gòu)的屈曲、過(guò)屈曲、振動(dòng)做了研究，給出熱、電載荷對(duì)壓電-功能梯度板的屈曲、振動(dòng)的影響。不少科研人員還研究了熱、機(jī)、電多種載荷耦合作用下壓電層合智能結(jié)構(gòu)的靜動(dòng)態(tài)響應(yīng)［11～13］。其中，Yang 等［13］基于高階剪切理論，用以微分求積和Galerkin 法為基礎(chǔ)的半解析方法，研究了各種邊界條件下，功能梯度層合板的非線性振動(dòng)問(wèn)題。作者同時(shí)考慮了溫度，壓電層施加電壓驅(qū)動(dòng)，但是沒(méi)有考慮非均勻升溫，指出在大變形下固有頻率是非?！懊舾械摹薄Ｓ跐龋?4］采用Airy 應(yīng)力函數(shù)法推導(dǎo)出了壓電功能梯度材料懸臂梁二維平面問(wèn)題解析解。丁皓江等［15，16］系統(tǒng)地研究了橫觀各向同性壓電板、殼的三維問(wèn)題，給出了相應(yīng)的解析通解。LI［17］等將功能梯度Timoshenko梁?jiǎn)栴}歸結(jié)為含有中面位移、橫截面轉(zhuǎn)角、橫截面內(nèi)力以及變形后的軸線弧長(zhǎng)共7個(gè)基本未知函數(shù)的非線性常微分方程邊值問(wèn)題，采用打靶法數(shù)值求解所得強(qiáng)非線性邊值問(wèn)題，研究了升溫時(shí)的屈曲和過(guò)屈曲響應(yīng)解，給出了梁的過(guò)屈曲平衡路徑，詳細(xì)分析了熱載荷和材料梯度性質(zhì)對(duì)變形及過(guò)屈曲的影響。

但作者未見(jiàn)關(guān)于壓電層合梁和壓電功能梯度材料梁的分析中，對(duì)熱-電載荷共同作用下FGM-Euler梁屈曲前的自由振動(dòng)問(wèn)題，給出半解析數(shù)值解的報(bào)道。本工作擬在文獻(xiàn)［18，19］的基礎(chǔ)上，對(duì)上下表面粘貼壓電層的FGM 層合梁在橫向非均勻升溫以及對(duì)壓電層施加驅(qū)動(dòng)電壓載荷下的幾何非線性大變形進(jìn)行定量分析?；谳S向可伸長(zhǎng)梁的精確幾何非線性理論，建立壓電-功能梯度Euler 梁在熱-電-機(jī)械載荷作用下的幾何非線性動(dòng)力學(xué)控制方程，將問(wèn)題歸結(jié)為含有中面位移、橫截面轉(zhuǎn)角、橫截面內(nèi)力以及變形后的軸線弧長(zhǎng)等14個(gè)基本未知函數(shù)的非線性常微分方程邊值問(wèn)題。把問(wèn)題線性化以后，得到兩端固定和簡(jiǎn)支的壓電-FGM 層合梁在橫向均勻電場(chǎng)和非均勻升溫作用下線性振動(dòng)的半解析數(shù)值解。

1 問(wèn)題的數(shù)學(xué)模型

考慮壓電材料和功能梯度材料組成的層合梁。其上、下表面層為壓電材料，厚度均為t，中間層為功能梯度材料，厚度為h。梁的長(zhǎng)度l，寬度為b，總厚度為h0。為研究方便，只考慮逆壓電效應(yīng)，并且認(rèn)為壓電層與功能梯度層之間牢固粘貼。這里，研究梁在溫度載荷和電壓作用下的靜動(dòng)態(tài)響應(yīng)，分析通過(guò)對(duì)壓電層的電壓驅(qū)動(dòng)實(shí)現(xiàn)對(duì)梁的熱靜態(tài)變形，振動(dòng)固有頻率的調(diào)整作用。在問(wèn)題的數(shù)學(xué)模型建立中考慮梁的軸向伸長(zhǎng)。

1.1 非線性幾何方程

精確考慮梁的軸線伸長(zhǎng)，可得梁變形后軸線的幾何方程［17，18］:

其中，x 為變形前梁軸線上任意一點(diǎn)的物質(zhì)坐標(biāo)，且x ∈［0，l］;u0(x)和w0(x)分別為軸線上的點(diǎn)在x，z 方向的位移，θ(x)為梁變形后軸線切線與x 軸正向的夾角，Λ0為軸線伸長(zhǎng)率。假設(shè)梁變形過(guò)程中橫截面仍然保持為平面，可以求得橫截面上任意點(diǎn)處的正應(yīng)變?yōu)椋?7～19］:

其中θ 為橫截面法線與x 軸之間的夾角。

1.2 本構(gòu)方程

假設(shè)材料為線性熱彈性材料，并且梁的各層之間截面上的變形是連續(xù)的，由此可得壓電-FGM 層合梁中第j 層的物理方程［9，10］:

其中:Ej，Gj，αj和E'zj分別為第j 層的拉伸彈性模量、剪切彈性模量、熱膨脹系數(shù)和電場(chǎng)強(qiáng)度(j = 1，2，3 )，d31為壓電應(yīng)變常量，T 為升溫?？紤]外電壓只沿厚度方向作用，則壓電層的電場(chǎng)強(qiáng)度可以表示為［9，10］:

Vj為第j 層的電壓，t1= t3= t，t2=h。這里由于只有上下表層存在電壓，所以:V2= 0 。

設(shè)功能梯度材料中間層由陶瓷和金屬兩種材料組成，并且由上表面的純陶瓷連續(xù)變化為下表面的純金屬。FGM 的材料的物性參數(shù)(彈性模量E、泊松比ν、熱膨脹系數(shù)α 和熱傳導(dǎo)系數(shù)K 等)滿足下列關(guān)系［8，10］

其中Vm，Vc分別為金屬和陶瓷的體積分?jǐn)?shù)。這里假設(shè)陶瓷體積分?jǐn)?shù)沿厚度方向坐標(biāo)按冪函數(shù)變化，即:

其中n 為陶瓷材料的體積分?jǐn)?shù)指數(shù)(0 ≤n ＜∞)?？紤]到功能梯度材料仍屬于各向同性材料，所以彈性常數(shù)之間滿足關(guān)系G = E/［2(1 + ν)］。

考慮到壓電層為各向同性材料;極化方向只沿z 軸方向。將式(2)代入式(3)可得內(nèi)力表達(dá)式:

其中C1，C2和C3分別為剛度系數(shù);NT，NE分別為熱軸力和電軸力，MT，ME分別為熱彎矩和電彎矩，具體定義為:

1.3 熱傳導(dǎo)方程

考慮層合梁受橫向非均勻升溫作用。設(shè)上表面升溫為T(mén)U，下表面升溫為T(mén)L。梁內(nèi)的升溫場(chǎng)可由下列一維熱傳導(dǎo)方程確定:

由于壓電層為均勻材料，熱傳導(dǎo)系數(shù)為常數(shù)。因此，由式(16)可知上、下層的升溫為線性函數(shù)。中間層功能梯度材料的熱傳導(dǎo)系數(shù)為橫向坐標(biāo)z 的已知函數(shù)，具體可由式(14)給出。由上、下表面的邊界條件和界面出處的連續(xù)性條件可以確定各層的升溫場(chǎng)分別如下［10］:

其中T12和T32分別為功能梯度層上表面和下表面(或上、下界面)處的溫度。具體表達(dá)式為:

其中Ki( i=p，c，m )分別表示壓電、陶瓷和金屬材料的導(dǎo)熱系數(shù);c 為常數(shù)，其表達(dá)式為:

其中

1.4 平衡方程

采用軸線可伸長(zhǎng)幾何非線性理論，可得Euler梁在變形后構(gòu)形的平衡方程［18，19］:

其中:H，V 分別為橫截面內(nèi)的水平和橫向內(nèi)力分量;(qx，qz)為作用于梁上的分布力;qθ為分布力矩。由靜力等效關(guān)系知，變形后橫截面內(nèi)法向內(nèi)力N 可由H，V 表示為:

由方程(7，8)和(21)可以解得:

若只考慮自由振動(dòng)，其慣性力可以表示為:

其中:I0，I1分別為軸線伸長(zhǎng)后單位長(zhǎng)度的質(zhì)量分布和單位長(zhǎng)度的轉(zhuǎn)動(dòng)慣性矩。表達(dá)式為:

至此，式(1)，(20)，(22)就構(gòu)成了壓電-功能梯度Euler-Bernoulli 梁在熱-電-機(jī)載荷同時(shí)作用下的幾何非線性力學(xué)行為的控制方程。其中包含了軸線位移u0(x)，w0(x)，橫截面轉(zhuǎn)角θ(x)以及等效內(nèi)力H(x)，V(x)和彎矩M(x)等14個(gè)基本未知函數(shù)，它們都是物質(zhì)坐標(biāo)x 的函數(shù)。由式(22)，(23)可知，軸線伸長(zhǎng)率Λ0和橫截面轉(zhuǎn)角θ 也可以表示為上述基本未知函數(shù)的表達(dá)式。

2 無(wú)量綱控制方程

為便于計(jì)算和分析，將上述控制方程轉(zhuǎn)化為無(wú)量綱形式。為此，引入下列無(wú)量綱量:

將上述無(wú)量綱變換代入式(1)，(20)，(22)，得到壓電-FGM 層合梁無(wú)量綱的動(dòng)力學(xué)控制方程:

其中:

控制方程中的無(wú)量綱電軸力、電彎矩的表達(dá)式見(jiàn)式(28)，(29)。上下表層電場(chǎng)同向(即V1= － V2= V［7］)時(shí):

上下表層電場(chǎng)異向(即V1= V2= V［7］)時(shí):

考慮梁的約束是兩端為固定和簡(jiǎn)支，則相應(yīng)的無(wú)量綱邊界條件可表示為式(30)和(31):

1)fixed-fixed

2)pinned-pinned

假設(shè)梁處于平衡狀態(tài)，在控制方程式(26)中，假設(shè)令sinθ = θ，cosθ = 1，并略去所有非線性項(xiàng)，則可得線性化問(wèn)題:

從上面的方程中消去PH，PV，m，θ，得到只含基本位移U，W 的方程:

于是得到對(duì)稱(chēng)貼壓電層的功能梯度梁在未屈曲狀態(tài)的線性振動(dòng)控制方程式(33)?？梢钥闯鰴M向振動(dòng)和縱向振動(dòng)是非耦合的，兩式可以獨(dú)立求解。設(shè)梁的自由振動(dòng)響應(yīng)為:

其中，ω1，ω2分別為縱向和橫向振動(dòng)的無(wú)量綱固有頻率，φ1，φ2為常數(shù)。將式(34a，b)帶入控制方程式(33)得到下列線性邊值問(wèn)題:

pinned-pinned:

fixed-fixed:

3 結(jié)果與討論

為了滿足u 的邊界條件，假設(shè)u = Asinnπξ，代入控制方程式(35a)，得到梁縱向振動(dòng)的各階固有頻率的解析表達(dá)式:

縱向振動(dòng)的一階固有頻率為:

為驗(yàn)證結(jié)果的正確性，將貼壓電層的功能梯度材料梁退化到各向同性材料(無(wú)壓電，功能梯度參數(shù))時(shí)，n = 0，F(xiàn)2= 0，F(xiàn)1= 12δ2，此時(shí)ω(n)1= nπ與文獻(xiàn)［24］中的結(jié)果相同。

當(dāng)梁的兩端簡(jiǎn)支pinned-pinned 時(shí)，w = Bsinnπξ滿足相應(yīng)的邊界條件，代入方程式(35b)，由B ≠0可得特征方程:

從而求得系統(tǒng)橫向振動(dòng)頻率

其中:μTE= F2(mT+ mE)－ (PT+ PE)。

由梁的熱彈性穩(wěn)定性分析可知，μ(n)= (nπ)2是對(duì)應(yīng)于pinned-pinned 梁的熱屈曲各階失穩(wěn)模態(tài)的載荷特征值(無(wú)量綱升溫)。n = 1 時(shí)得到臨界升溫μcr= μ(1)= π2，于是式(38)可以表示為:

如果退化到均勻(陶瓷)材料，則m1= m2= 0;PE= 0;mE= 0;n = 0 ;F1=12δ2，F(xiàn)2= 0 ，F(xiàn)0=1/12δ2，PT= 12δ2αcT，于 是 得 到:=則與橫向振動(dòng)頻率的表達(dá)式其中μ = 12δ2αcT，r2= 1/12δ2比較，可以看出在均勻升溫下兩者的表達(dá)式是完全一樣的。

對(duì)于fixed-fixed 梁，選取滿足相應(yīng)橫向邊界條件的近似的一階振型函數(shù):

利用Galerkin 法［5，13］，可以得到梁的一階近似固有頻率

對(duì)于其他低階固有頻率則考慮打靶法［17～19］求解，計(jì)算中，考慮梁的中間層為陶瓷氧化鋯(Zirconia/ZrO2)和金屬鋁(Aluminum/Al)兩相材料制成的功能梯度材料(記為ZrO2/Al)，上、下表層壓電材料為G-1195N。材料物性參數(shù)見(jiàn)表1［9，10］。

表1 壓電-功能梯度層合梁的物性參數(shù)Table 1 Material property parameters of the piezoelectric FGM beams

為了驗(yàn)證上述理論分析和計(jì)算程序的正確性，在將壓電-功能梯度梁退化為均勻各向同性梁的情況下(n = 0 ，η = 0 )的所得數(shù)值結(jié)果與已有結(jié)果進(jìn)行了比較。在均勻升溫情況下(Tr= 1 )，表2 中給出由本工作所計(jì)算獲得的無(wú)量綱臨界屈曲溫度τcr= 12δ2αcTcr與文獻(xiàn)［17］中結(jié)果的比較。可以看出，兩者吻合較好。

對(duì)于控制方程(35b)和邊界條件(35d)，利用打靶法把邊值轉(zhuǎn)化為初值問(wèn)題求解，退化到無(wú)壓電，均勻材料，均勻升溫情況下，求得在臨界溫度附近的其他低階固有頻率:

表2 具有不同長(zhǎng)細(xì)比梁的無(wú)量綱臨界升溫值Table 2 Critical non-dimensional thermal of beams with aspect radio δ

由頻率的表達(dá)式(39)，(42)知，對(duì)于細(xì)長(zhǎng)梁，截面轉(zhuǎn)動(dòng)慣量(含在F0里)對(duì)橫向振動(dòng)頻率可以忽略，在不考慮升溫、壓電情況下，表3 給出了細(xì)長(zhǎng)比δ 為不同值時(shí)，考慮轉(zhuǎn)動(dòng)慣性力時(shí)的一階固有頻率與忽略轉(zhuǎn)動(dòng)慣性力時(shí)的一階頻率的比值，由此可見(jiàn)轉(zhuǎn)動(dòng)慣性力對(duì)頻率的影響程度。

表3 轉(zhuǎn)動(dòng)慣性力對(duì)振動(dòng)頻率值比的影響(PT，PE，MT，ME = 0 )Table 3 Effects of inertia forces to natural frequencies

表3 轉(zhuǎn)動(dòng)慣性力對(duì)振動(dòng)頻率值比的影響(PT，PE，MT，ME = 0 )Table 3 Effects of inertia forces to natural frequencies

δ Pinned-pinned Fixed-fixed 5 0.9839 0.9945 10 0.9959 0.9986 15 0.9981 0.9993 20 0.9989 0.9996 30 0.9995 0.9998 40 0.9997 0.9999 50 0.9998 0.9999 100 0.9999 0.9999

表4 ～表7 給出了固支條件下，不同載荷參數(shù)、不同材料參數(shù)情況下，F(xiàn)GM 層合梁在臨界載荷附近的其他低階固有頻率。

4 結(jié)論

(1)梁縱向振動(dòng)頻率只與材料的性質(zhì)有關(guān)，與電壓升溫都沒(méi)有關(guān)系，這是因?yàn)闇囟葓?chǎng)和電場(chǎng)都是橫向的，對(duì)縱向振動(dòng)沒(méi)有影響。

(2)當(dāng)上下壓電層電壓同向時(shí)，升壓都會(huì)降低梁橫向振動(dòng)的頻率;當(dāng)上下壓電層電壓異向時(shí)，升壓都會(huì)增大梁橫向振動(dòng)頻率。

表4 在參數(shù)η = 0，Tr = 1，δ = 30 一定的情況下其他低階固有頻率Table 4 The other natural frequencies when parameters η = 0，Tr = 1，δ = 30 are given

表5 在參數(shù)λ = 0，Tr = 1，δ = 30 一定的情況下其他低階固有頻率Table 5 The other natural frequencies when parameters λ = 0，Tr = 1，δ = 30 are given

表6 在參數(shù)n = 0，Tr = 5，δ = 30 一定的情況下其他低階固有頻率Table 6 The other natural frequencies when parameters n = 0，Tr = 5，δ = 30 are given

表7 在參數(shù)λ = 0，Tr = 5，δ = 30 一定的情況下其他低階固有頻率Table 7 The other natural frequencies when parameters λ = 0，Tr = 5，δ = 30 are given

(3)升溫只會(huì)降低梁的固有頻率;溫度場(chǎng)非均勻程度也會(huì)影響結(jié)構(gòu)的固有頻率。

(4)通過(guò)對(duì)上下壓電層的反向驅(qū)動(dòng)(Φ ＜0 )，可以減小甚至抵消掉升溫引起的自振頻率的變化，實(shí)現(xiàn)對(duì)結(jié)構(gòu)的固有頻率的調(diào)節(jié)。

［1］CEAWLEY E F，de LUIS J. Use of piezoelectric actuators as elements of intelligent structures［J］. AIAA Journal，1987，25(10):1373 －1385.

［2］ZHOU Y H，WANG J Z. Vibration control of piezoelectric beam-type-plates with geometrical nonlinear deformation［J］. Int J Non-Linear Mech，2004，39:909 －920.

［3］XIA X K，SHEN H S. Vibration of postbuckled FGM hybrid laminated plates in thermal environment［J］. Engineering Structures，2008，30(9):2420 －2435.

［4］丁麗霞，劉瑋. 壓電元件驅(qū)動(dòng)的功能梯度彈性薄板的屈曲［J］. 功能材料，2006，37(8):1229 －1238.(DING L X，LIU W. The buckling analysis of functionally gradient thin elastic plate bonded piezoelectric patches［J］.Journal of Functional Materials，2006，37(8):1229 －1238.)

［5］DAI K Y，LIU G R ，HAN X，et al. Thermo-mechanical analysis of functionally graded material (FGM)plates using element-free Galerkin method［J］. Computers and Structures，2005，83:1487 －1502.

［6］林啟榮，劉正興，王宗利.電場(chǎng)作用下壓電層合梁的分析［J］.應(yīng)用數(shù)學(xué)和力學(xué)，2001，22 (9):969 －975.(LIN Q R，LIU Z X，WANG Z G. Analysis of beams with piezoelectric actuators［J］. Applied Mathematics and Mechanics，2001，22(9):969 －975.)

［7］FRIDMAN Y，ABRAMOVICH H. Enhanced structural behavior of flexible laminated composite beams［J］.Composite Structures，2008，82 :140 －154.

［8］SHEN H S. Postbuckling of FGM plates with piezoelectric actuators under thermo-electro-mechanical loadings［J］. Int J of Solids and Structures，2005，42:6101 －6121.

［9］SHEN H S. Post-buckling of shear deformable laminated with piezoelectric actuators under complex loading conditions［J］. Int J of Solids and Structures，2001，38:7703 －7721.

［10］HUANG X L，SHEN H S. Vibration and dynamic response of functionally graded plates with piezoelectric actuators in thermal environments［J］. Journal of Sound and Vibration，2006，289(1/2):25 －53.

［11］ANNIGERI A R，GANESAN N，SWARNAMANI S. Free vibration behaviour of multiphase and layered magneto-electro-elastic beam［J］. Journal of Sound and Vibration，2007，299 :44 －63.

［12］RAJESH K，BHANGALE，GANESAN N. Free vibration studies of simply supported non-homogeneous functionally graded magneto-electro-elastic finite cylindrical shells［J］.Journal of Sound and Vibration ，2005，288:412 －422.

［13］YANG J，KITIPORNCHAI S，LIEW K M. Large amplitude vibration of thermo-electro-mechanically stressed FGM laminated plates［J］. Comput Methods Appl Mech Engrg，2003，192:3861 －3885.

［14］于濤，仲政.功能梯度壓電懸臂梁的彎曲分析［J］.中國(guó)科學(xué):物理學(xué) 力學(xué) 天文學(xué)，2006，36(5):518 －529.(YU T，ZHONG Z. Bending analysis of a functionally graded piezoelectric cantilever beam［J］. Scientia Sinica:Physica，Mechanica ＆ Astronomica，2006，36(5):518 －529.)

［15］丁皓江，徐榮橋，國(guó)風(fēng)林. 軸對(duì)稱(chēng)層合橫觀各向同性壓電圓板的精確解［J］. 中國(guó)科學(xué):技術(shù)科學(xué)，1999，29(3):214 －221.(DING H J，XU R Q，GUO F L.Exact solution of axisymmetric laminated transversely isotropic piezoelectric circular plate［J］. Scientia Sinica:Technologica，1999，29 (3):214 －221.)

［16］丁皓江，陳偉球，徐榮橋. 壓電板殼自由振動(dòng)的三維精確分析［J］.力學(xué)季刊，2001，22 (1):1 －9.(DING H J，CHEN W Q，XU R Q. Three-dimensional exact analyses of free vibrations for piezoelectric plates and shell［J］. Chinese Quarterly of Mechanics，2001，22(1):1 －9.)

［17］LI S R，ZHANG J H，ZHAO Y G. Thermal post-buckling of functionally graded material timoshenko beams［J］. Applied Mathematics and Mechanics (English Edition)，2006，27(6):709 －715.

［18］LI S R，SU H D，CHENG C J. Free vibration of functionally graded material beams with surface-bonded piezoelectric layers in thermal environment［J］. Applied Mathematics and Mechanics(English Edition)，2009，30(8):1 －14.

［19］蘇厚德，李世榮，高穎. 粘貼壓電層功能梯度材料Timoshenko 梁的熱過(guò)屈曲分析［J］. 計(jì)算力學(xué)學(xué)報(bào)，2010，27(6):1067 －1072.(SU H D，LI S R，GAO Y，Thermal post-buckling of functionally graded material timoshenko beams with surface-bonded piezoelectric layers［J］. Chinese Journal of Computational Mechanics，2010，27 (6 ):1067 －1072. )