杜超凡 鄭燕龍 周曉婷 章定國

*(揚(yáng)州大學(xué)建筑科學(xué)與工程學(xué)院，江蘇揚(yáng)州 225000)

?(南京理工大學(xué)理學(xué)院，南京 210094)

功能梯度材料 (functionally graded materials，F(xiàn)GM) 是由兩種或兩種以上性能各異的材料組合而成，通過改變其體積分?jǐn)?shù)來改變結(jié)構(gòu)沿某一方向的性能，以實(shí)現(xiàn)結(jié)構(gòu)不同位置的功能需求[1]。由于FGM在材料界面處的力學(xué)性能過渡平滑，能消除不必要的界面應(yīng)力，為解決熱障問題提供了有力的幫助，因此其性能優(yōu)于由兩種離散材料粘結(jié)在一起的復(fù)合增強(qiáng)材料，被研究人員廣泛關(guān)注。這種固有的材料性能不均勻性給力學(xué)分析帶來了很大的困難，以往針對(duì)各向同性、均勻材料引入和發(fā)展的力學(xué)概念和建模理論等已不再適用，需要進(jìn)行探索和創(chuàng)新。由FGM制成的旋轉(zhuǎn)葉片，如太陽能帆板、風(fēng)力發(fā)電機(jī)葉片及航空發(fā)動(dòng)機(jī)葉片等，其動(dòng)力學(xué)特性必然有別于傳統(tǒng)均質(zhì)材料。因此研究該類結(jié)構(gòu)在不同材料組分和梯度分布規(guī)律下的動(dòng)力學(xué)特性，將為工程應(yīng)用設(shè)計(jì)出理想的功能梯度材料提供理論指導(dǎo)，具有重要的理論意義和廣闊的應(yīng)用前景。

旋轉(zhuǎn)葉片是典型的剛?cè)狁詈戏蔷€性系統(tǒng)，越來越多的學(xué)者關(guān)注柔性葉片的大范圍運(yùn)動(dòng)與自身彈性變形的剛?cè)狁詈蠙C(jī)理以及柔性體變形場(chǎng)的離散方法，對(duì)這些問題的研究已成為柔性多體系統(tǒng)動(dòng)力學(xué)領(lǐng)域的熱點(diǎn)和難點(diǎn)。自從Kane 等[2]揭示“動(dòng)力剛化”現(xiàn)象以來，國內(nèi)外許多學(xué)者均對(duì)板的剛?cè)狁詈蠙C(jī)理進(jìn)行了多方位研究。劉錦陽等[3]從連續(xù)介質(zhì)力學(xué)出發(fā)，基于Jourdan 速度變分原理，采用有限元法和假設(shè)模態(tài)法對(duì)均質(zhì)材料旋轉(zhuǎn)矩形薄板動(dòng)力學(xué)響應(yīng)進(jìn)行了仿真，并對(duì)兩種方法的結(jié)果進(jìn)行了對(duì)比。Yoo 等[4]基于Kirchhoff假設(shè)，采用Rayleigh–Ritz 方法描述三個(gè)維度的變形，考慮縱向和橫向變形耦合及運(yùn)動(dòng)引起的動(dòng)力剛化效應(yīng)，研究了旋轉(zhuǎn)復(fù)合材料懸臂板的不同纖維角對(duì)固有頻率和模態(tài)的影響。Niu 等[5]運(yùn)用哈密頓原理結(jié)合一階剪切變形理論和馮·卡門幾何非線性研究了受橫向激勵(lì)作用下的簡(jiǎn)支功能梯度石墨烯復(fù)合材料圓柱板的非線性振動(dòng)問題。Javani 等[6]基于Mindlin 板理論，采用廣義微分求積法(generalized differential quadrature method，GDQM) 研究了功能梯度環(huán)形薄板的熱致振動(dòng)問題。Alireza 等[7]通過Kirchhoff假設(shè)和一階剪切變形理論確立非線性方程和邊界條件，運(yùn)用GDQM 和Newton–Raphson迭代法對(duì)方程離散和求解，基于Crank–Nicolson 方案獲得溫度場(chǎng)的函數(shù)，研究了FGM 薄板遭遇溫度突兀冷卻的非線性熱彈性瞬時(shí)響應(yīng)。Li 等[8]和黎亮等[9]基于假設(shè)模態(tài)法，考慮不同功能梯度系數(shù)對(duì)動(dòng)力方程的影響，先后研究了FGM 矩形薄板自由振動(dòng)響應(yīng)及彎曲和扭轉(zhuǎn)之間的耦合模態(tài)和頻率轉(zhuǎn)向現(xiàn)象。劉璟澤等[10]采用Kirchhoff–Mindlin 理論假設(shè)剪切應(yīng)變場(chǎng)規(guī)避剪切閉鎖，使用有限元三角形單元進(jìn)行變形場(chǎng)描述，對(duì)曲線加筋Kirchhoff–Mindlin 板自由振動(dòng)進(jìn)行了分析。Moita 等[11]使用有限元法分析了復(fù)合板FGM 層、壓電層和具有黏彈性的核心夾層之間主動(dòng)-被動(dòng)阻尼[12]的振動(dòng)問題，F(xiàn)GM 層和壓電層皆用經(jīng)典薄板理論建模，核心夾層使用三階剪切變形理論建模。Mantari 等[13]基于只包含四個(gè)未知量的新型一階剪切變形理論，通過哈密頓原理推導(dǎo)控制方程，對(duì)FGM 板進(jìn)行了自由振動(dòng)分析。Zhang等[14]采用攝動(dòng)法，將壓電復(fù)合材料矩形薄板在橫向和面內(nèi)激勵(lì)下的復(fù)雜偏微分非線性運(yùn)動(dòng)控制方程轉(zhuǎn)化為等效的可解非線性方程，進(jìn)而分析了矩形薄板的運(yùn)動(dòng)響應(yīng)。Subodh 等[15]基于動(dòng)態(tài)剛度法對(duì)功能梯度矩形薄板的自由振動(dòng)行為進(jìn)行了分析。吳根勇等[16]基于經(jīng)典薄板理論，考慮二次耦合變形量，采用有限元法對(duì)復(fù)合材料板進(jìn)行離散，建立了作大范圍運(yùn)動(dòng)復(fù)合材料板的動(dòng)力學(xué)方程，研究了材料鋪層角度對(duì)復(fù)合材料變形的影響并和同性材料進(jìn)行了差異對(duì)比。Li 等[17]和Akhras 等[18]均使用B 樣條法對(duì)變形場(chǎng)離散，分別對(duì)壓電層復(fù)合板進(jìn)行了方法優(yōu)越性分析和穩(wěn)定性分析。楊興等[19]從一階剪切變形理論出發(fā)，使用有限元法對(duì)變形場(chǎng)進(jìn)行離散，研究了作大范圍運(yùn)動(dòng)功能梯度材料厚板的響應(yīng)，并同經(jīng)典薄板理論進(jìn)行比較，驗(yàn)證了經(jīng)典薄板理論的不足之處。然而，采用新的變形場(chǎng)離散方法，對(duì)考慮板剪切變形的旋轉(zhuǎn)FGM 板動(dòng)力學(xué)問題的研究比較少見。

本文采用無網(wǎng)格徑向基點(diǎn)插值法[20-21]描述板的變形場(chǎng)，考慮材料非均勻性以及橫向彎曲引起的縱向縮短，即非線性耦合項(xiàng)，利用浮動(dòng)坐標(biāo)系法[22]，基于一階剪切變形理論，運(yùn)用第二類拉格朗日方程建立了作旋轉(zhuǎn)運(yùn)動(dòng)FGM 矩形板的動(dòng)力學(xué)方程。研究了作旋轉(zhuǎn)運(yùn)動(dòng)功能梯度材料矩形板的動(dòng)力學(xué)行為隨功能梯度指數(shù)、轉(zhuǎn)速、板縱橫比的變化情況。同時(shí)對(duì)勻速轉(zhuǎn)動(dòng)的中心剛體與FGM 矩形板系統(tǒng)自由振動(dòng)固有頻率的頻率轉(zhuǎn)向問題進(jìn)行了研究。

1 動(dòng)力學(xué)建模

1.1 功能梯度材料

功能梯度材料通常是由兩種材料混合而成，一種是陶瓷，一種是金屬。對(duì)于FGM 矩形彈性板而言，假定其板上表面為陶瓷，下表面為金屬，兩表面的中間部分由這兩種材料組合物構(gòu)成。本文FGM 矩形板的材料表示如下：彈性模量E(z)，材料的密度ρ(z)，泊松比μ(z)，可由下式表示[23]

式中P代表 FGM 矩形板某一位置的物理特性E(z)，ρ(z)，μ(z)，Pc和Pm分別代表FGM 矩形板上表面z=h/2 和下表面z=-h/2 陶瓷材料和金屬材料的物理特性，h為板厚度，V代表陶瓷材料的體積分?jǐn)?shù)，N(N≥0) 為FGM 從上表面過渡到下表面的梯度指數(shù)。

1.2 運(yùn)動(dòng)學(xué)描述

圖1 為在空間中作大范圍運(yùn)動(dòng)的FGM 矩形板模型。其中O-XY Z為慣性坐標(biāo)系，o-xyz為連體坐標(biāo)系，且連體坐標(biāo)系所在平面x-o-y與FGM 板未變形前中面重合，e1，e2，e3分別為連體坐標(biāo)三個(gè)坐標(biāo)軸的單位矢量。FGM 板長為L，寬為H，厚度為h，密度為ρ(z)，彈性模量為E(z)，泊松比為μ(z)=μc=μm=μ。

圖1 作大范圍運(yùn)動(dòng)FGM 矩形板

R為中心剛體半徑，r0為連體坐標(biāo)系原點(diǎn)在慣性坐標(biāo)系中的矢徑，ρ0為點(diǎn)P0在連體坐標(biāo)系中的位置?，F(xiàn)假定P0點(diǎn)為矩形板變形前板內(nèi)任意一點(diǎn)，變形后為P，變形矢量為u，在連體坐標(biāo)系下的分量為(u1，u2，u3)。其中縱向位移u1，u2可表示為

式中，VO和ωA分別為連體坐標(biāo)系相對(duì)于慣性坐標(biāo)O-XY Z的速度和角速度矢量，VP A為P點(diǎn)相對(duì)連體坐標(biāo)系的速度矢量。則P點(diǎn)相對(duì)于慣性坐標(biāo)系O-XY Z的速度矢量為

其中κ為剪切修正系數(shù)，取值為κ=5/6。

1.3 系統(tǒng)的動(dòng)能和勢(shì)能

FGM 矩形板的動(dòng)能為

將式(5) 和式(6) 代入式(9) 得

1.4 FGM 矩形板變形場(chǎng)的離散

常用的矩形板變形場(chǎng)離散方法有假設(shè)模態(tài)法和有限元法，而本文采用無網(wǎng)格法中的徑向基點(diǎn)插值法(radial point interpolation method， RPIM)。對(duì)于問題域Ω中的連續(xù)場(chǎng)函數(shù)u(x)，可以用在該問題域中相關(guān)點(diǎn)的基函數(shù)表示如下

式中dc為計(jì)算點(diǎn)x同其局部支持域中的節(jié)點(diǎn)的間距的均值。

矩形板問題域中任意一點(diǎn)的變形可以表示為

式中下標(biāo)“，” 表示對(duì)坐標(biāo)求偏導(dǎo)數(shù)。

1.5 FGM 矩形板剛?cè)狁詈蟿?dòng)力學(xué)方程

得到FGM 板的剛?cè)狁詈蟿?dòng)力學(xué)方程為

式中各分塊矩陣為

式中，a01，a02，a03為基點(diǎn)加速度在連體坐標(biāo)系下的分量，表達(dá)式為

各分塊式中雙下劃線項(xiàng)為考慮二次耦合變形量帶來的附加動(dòng)力剛度項(xiàng)。傳統(tǒng)的零次近似模型建模時(shí)忽略了這些項(xiàng)，即這些項(xiàng)都為零。單下劃線為考慮橫向剪切變形所推導(dǎo)出的剪切變形項(xiàng)，若不考慮該項(xiàng)，則動(dòng)力學(xué)方程退化為經(jīng)典薄板理論。

2 數(shù)值仿真

2.1 功能梯度板“動(dòng)力剛化” 研究

本節(jié)仿真一作大范圍已知運(yùn)動(dòng)的懸臂FGM 矩形板， hub-FGM 板系統(tǒng)的材料參數(shù)： 長L=1.828 8 m，寬H=1.219 2 m，厚度h=0.02 m，ρc=3000 kg/m3，ρm=2707 kg/m3，Ec=151 GPa，Em=70 GPa，泊松比μ=0.3，中心剛體半徑R=0。矩形板節(jié)點(diǎn)劃分沿x和y軸分別以16 和8 個(gè)節(jié)點(diǎn)均勻分布，積分網(wǎng)格沿x和y軸分別均勻劃分8 格和6 格，假設(shè)模態(tài)法(assumed mode method，AMM) 和三角形單元有限元法(triangle finite element method，TRFEM) 分別取5×7 階模態(tài)和2×16×8 個(gè)單元。給定角速度規(guī)律為：ω1=ω3= ˙ω1= ˙ω3= 0，ω2=ω，˙ω2= ˙ω。其中

式中，時(shí)間T=30 s。

若只考慮FGM 板沿厚度方向上的橫向振動(dòng)，系統(tǒng)的動(dòng)力學(xué)方程可簡(jiǎn)化為

圖2 和圖3 是角速度分別取Ω=3 Hz 和10 Hz，板功能梯度系數(shù)N分別取1 和2 的零次近似模型和一次近似模型計(jì)算得到的板外側(cè)角點(diǎn)的橫向變形圖。從圖2 中可以看出，本文模型的計(jì)算結(jié)果與文獻(xiàn)[9] 基本吻合，但是在最大變形處，計(jì)算結(jié)果存在微小差異。這是因?yàn)槲墨I(xiàn)[9] 是基于經(jīng)典薄板理論，而本文是基于考慮剪切變形的Mindlin 板理論，因而系統(tǒng)更顯柔性，所以最大變形處偏大。比較圖3(a)和圖3(b) 可以看出，Ω=3 Hz 低速運(yùn)動(dòng)時(shí)，不考慮附加剛度項(xiàng)的零次近似模型計(jì)算結(jié)果收斂，當(dāng)運(yùn)動(dòng)為高速運(yùn)動(dòng)Ω= 10 Hz 時(shí)，零次近似模型的計(jì)算結(jié)果則發(fā)散，而一次近似模型的結(jié)果依然收斂。說明未考慮附加剛度項(xiàng)的零次近似模型只適用于低轉(zhuǎn)速工況，而一次近似模型既適用于低速，又適用于高轉(zhuǎn)速工況。

圖2 N =1，Ω =3 Hz 時(shí)矩形板外側(cè)角點(diǎn)橫向變形

圖3 N =2 時(shí)矩形板外側(cè)角點(diǎn)橫向變形

圖4 和圖5 是N= 2，角速度Ω= 3 Hz 和10 Hz 時(shí)一次近似模型得到板外側(cè)角點(diǎn)橫向變形速度。從兩圖中可以看出三種離散方法的仿真結(jié)果基本一致。進(jìn)一步說明了本文所建模型的正確性。

圖4 N =2，Ω =3 Hz 時(shí)板外側(cè)角點(diǎn)橫向變形速度

圖5 N =2，Ω =10 Hz 時(shí)板外側(cè)角點(diǎn)橫向變形速度

圖6 為Ω= 10 Hz，h= 0.002 m，Ec=15.1 GPa，N= 1，其余參數(shù)不變的FGM 矩形板外側(cè)角點(diǎn)橫向變形圖。從圖中可以看出，TR-FEM和RPIM 的結(jié)果在前10 s 基本一致，但隨著時(shí)間的增大，TR-FEM 結(jié)果發(fā)散。而RPIM 的結(jié)果仍然收斂，說明了RPIM 有更好的計(jì)算穩(wěn)定性，在同等計(jì)算條件下比TR-FEM 更有計(jì)算優(yōu)勢(shì)。

圖6 矩形板外側(cè)角點(diǎn)橫向變形圖

圖7 分別給出了在Ω= 3 Hz 和Ω= 10 Hz 時(shí)板外側(cè)角點(diǎn)的橫向變形隨功能梯度系數(shù)的變化情況。從圖中可以看出，當(dāng)N= 0 時(shí)，功能梯度板退化為均質(zhì)板，板的最大變形最小，隨著N的增大，板外側(cè)角點(diǎn)變形也增大，說明功能梯度系數(shù)能影響板的柔性，功能梯度系數(shù)越大，板的柔性越大。

圖7 不同功能梯度系數(shù)下矩形板外側(cè)角點(diǎn)的橫向變形

2.2 功能梯度板的固有頻率研究

假設(shè)中心剛體半徑為R，中心剛體以角速度Ω繞y軸勻速轉(zhuǎn)動(dòng)，則在連體坐標(biāo)系下，原點(diǎn)O的速度矢量為：(v1，v2，v3) = (0，0，-RΩ)，(ω1，ω2，ω3) =(0，Ω，0)。

若只考慮矩形板的橫向振動(dòng)，則基于Mindlin 板理論推導(dǎo)的作定軸轉(zhuǎn)動(dòng)矩形板的自由振動(dòng)方程為

考慮到矩形板幾何尺寸對(duì)固有頻率的影響，定義如下無量綱變量：τ=t/T，χ=x/L，ζ=y/H，Z=q3/L，δ=L/H，η=h/L，σ=R/L，λ=Ec/Em，?=ρc/ρm，?=ωT，γ=ΩT。

將上述無量綱變量代入式(62) 得

為進(jìn)一步驗(yàn)證本文所建模型的正確性，將功能梯度板退化為勻質(zhì)薄板，同文獻(xiàn)[25] 中的計(jì)算結(jié)果進(jìn)行比較。表1 為δ=1，η=0.01，σ=0，中心剛體懸臂板無量綱角速度γ=1 的前五階無量綱固有頻率。由表1 可知，本文結(jié)果同文獻(xiàn)[25]基本吻合，且由于考慮剪切變形，系統(tǒng)更偏柔性，固有頻率比文獻(xiàn)中的結(jié)果都偏小。表2 為δ=1，η=0.01，σ=0，N=2，不同無量綱角速度γ下中心剛體FGM 懸臂板前三階無量綱固有頻率。由該表知，隨著無量綱角速度的增大，三種方法的無量綱固有頻率皆增大，且基于一階剪切變形理論的RPIM 結(jié)果總是比基于經(jīng)典薄板理論的假設(shè)模態(tài)法的數(shù)值小，說明忽略剪切變形的經(jīng)典薄板理論總是高估了系統(tǒng)的固有頻率。表3 為δ= 1，η= 0.01，σ= 0，γ= 10，N取不同值時(shí)，中心剛體FGM 懸臂板前五階無量綱固有頻率。結(jié)果表明，隨著N的增大，系統(tǒng)無量綱固有頻率減小，進(jìn)一步說明功能梯度系數(shù)越大，系統(tǒng)的柔性越大。

表1 作定軸轉(zhuǎn)動(dòng)中心剛體懸臂板的前五階無量綱固有頻率(δ =1，η =1，σ =0，γ =1)

表2 作定軸轉(zhuǎn)動(dòng)中心剛體FGM 懸臂板的前三階無量綱固有頻率(δ =1，η =0.01，σ =0，N =2)

表3 作定軸轉(zhuǎn)動(dòng)中心剛體FGM 懸臂板的前五階無量綱固有頻率(δ =1，η =0.01，σ =0，γ =10)

圖8 所示為δ=5，η=0.01，σ=1 時(shí)，均質(zhì)懸臂板前五階無量綱固有頻率隨無量綱角速度的變化圖。從圖中可以看出，各階無量綱固有頻率皆隨無量綱角速度的增大而增大，而且在第二階和第三階無量綱固有頻率有頻率轉(zhuǎn)向現(xiàn)象。從放大圖8(b) 可以看出，第二階和第三階固有頻率沒有發(fā)生交叉。

圖8 均質(zhì)懸臂板前五階無量綱固有頻率(δ =5，η =0.01，σ =1)

圖9 為δ=5，η=0.01，σ=1，N=2 時(shí)，F(xiàn)GM懸臂板前八階無量綱固有頻率隨無量綱角速度變化圖。從圖9(a)中能觀察到隨著無量綱角速度的增大，多階無量綱固有頻率發(fā)生了顯著的頻率轉(zhuǎn)向現(xiàn)象，而從頻率轉(zhuǎn)向區(qū)域的局部放大圖9(b)～9(d)中可知，發(fā)生頻率轉(zhuǎn)向的相鄰兩階頻率并沒有交叉。

圖9 FGM 懸臂板前八階無量綱固有頻率的變化情況(δ =5，η =0.01，σ =1，N =2)(續(xù))

圖10 為δ=1，η=0.01，σ=1，N=2 時(shí)，F(xiàn)GM懸臂板第六階到第八階無量綱固有頻率隨無量綱角速度變化圖和局部放大圖，從圖10(a)看到第七階固有頻率同第六階和第八階固有頻率發(fā)生了多次轉(zhuǎn)向現(xiàn)象，說明板的縱橫比對(duì)系統(tǒng)固有頻率的頻率轉(zhuǎn)向影響較大。

圖10 FGM 懸臂板第六階到第八階無量綱固有頻率的變化情況(δ =1，η =0.01，σ =1，N =2)

3 結(jié)論

(1)采用RPIM 描述柔性板的變形場(chǎng)，基于一階剪切變形理論，建立了旋轉(zhuǎn)FGM 矩形板的剛?cè)狁詈弦淮谓颇Ｐ?。該模型既適用于低轉(zhuǎn)速情況，又適用于高轉(zhuǎn)速情況，既適用于薄板問題，又適用于中厚板問題。在同等計(jì)算條件下，RPIM 相比TR-FEM 更具有計(jì)算優(yōu)越性。

(2) 隨著功能梯度指數(shù)的增大，旋轉(zhuǎn)FGM 矩形板的橫向彎曲變形變大，固有頻率減小，說明功能梯度指數(shù)的增大使系統(tǒng)的柔性增大，可通過改變功能梯度指數(shù)的取值范圍以達(dá)到控制變形及頻率的目的。

(3)基于經(jīng)典薄板理論的無量綱固有頻率總是大于基于一階剪切變形理論的無量綱固有頻率，說明忽略剪切變形的Kirchhoff假設(shè)使結(jié)構(gòu)更偏剛性，總是高估了系統(tǒng)的固有頻率。

(4)系統(tǒng)的無量綱固有頻率會(huì)隨著轉(zhuǎn)速和板縱橫比的變化而發(fā)生頻率轉(zhuǎn)向現(xiàn)象。不僅有單次轉(zhuǎn)向，還有多次轉(zhuǎn)向現(xiàn)象，且發(fā)生頻率轉(zhuǎn)向的相鄰兩階頻率并不相交。