作旋轉運動功能梯度材料矩形Mindlin 板的剛柔耦合動力學建模1)

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

*(揚州大學建筑科學與工程學院，江蘇揚州 225000)

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

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

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

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

1 動力學建模

1.1 功能梯度材料

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

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

1.2 運動學描述

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

圖1 作大范圍運動FGM 矩形板

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

式中，VO和ωA分別為連體坐標系相對于慣性坐標O-XY Z的速度和角速度矢量，VP A為P點相對連體坐標系的速度矢量。則P點相對于慣性坐標系O-XY Z的速度矢量為

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

1.3 系統(tǒng)的動能和勢能

FGM 矩形板的動能為

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

1.4 FGM 矩形板變形場的離散

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

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

矩形板問題域中任意一點的變形可以表示為

式中下標“，” 表示對坐標求偏導數(shù)。

1.5 FGM 矩形板剛柔耦合動力學方程

得到FGM 板的剛柔耦合動力學方程為

式中各分塊矩陣為

式中，a01，a02，a03為基點加速度在連體坐標系下的分量，表達式為

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

2 數(shù)值仿真

2.1 功能梯度板“動力剛化” 研究

本節(jié)仿真一作大范圍已知運動的懸臂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é)點劃分沿x和y軸分別以16 和8 個節(jié)點均勻分布，積分網格沿x和y軸分別均勻劃分8 格和6 格，假設模態(tài)法(assumed mode method，AMM) 和三角形單元有限元法(triangle finite element method，TRFEM) 分別取5×7 階模態(tài)和2×16×8 個單元。給定角速度規(guī)律為：ω1=ω3= ˙ω1= ˙ω3= 0，ω2=ω，˙ω2= ˙ω。其中

式中，時間T=30 s。

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

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

圖2 N =1，Ω =3 Hz 時矩形板外側角點橫向變形

圖3 N =2 時矩形板外側角點橫向變形

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

圖4 N =2，Ω =3 Hz 時板外側角點橫向變形速度

圖5 N =2，Ω =10 Hz 時板外側角點橫向變形速度

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

圖6 矩形板外側角點橫向變形圖

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

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

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

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

若只考慮矩形板的橫向振動，則基于Mindlin 板理論推導的作定軸轉動矩形板的自由振動方程為

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

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

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

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

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

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

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

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

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

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

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

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

3 結論

(1)采用RPIM 描述柔性板的變形場，基于一階剪切變形理論，建立了旋轉FGM 矩形板的剛柔耦合一次近似模型。該模型既適用于低轉速情況，又適用于高轉速情況，既適用于薄板問題，又適用于中厚板問題。在同等計算條件下，RPIM 相比TR-FEM 更具有計算優(yōu)越性。

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

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

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