郭琛琛， 劉 濤， 王青山， 秦 斌

(1. 中南大學(xué) 輕合金研究院， 長沙 410083； 2. 中南大學(xué) 高性能復(fù)雜制造國家重點實驗室， 長沙 410083；3. 中南大學(xué) 交通運輸工程學(xué)院， 長沙 410075)

復(fù)合材料層合板由于其較高的比強度和比剛度等力學(xué)特性，被廣泛應(yīng)用于土木工程、海洋工程和航空航天等諸多工業(yè)領(lǐng)域。這類構(gòu)件通常面臨各種極端惡劣的工況，并且可能受到形式各異的邊界載荷約束。因此，不同邊界條件下復(fù)合材料層合板的振動特性預(yù)測一直是研究人員關(guān)注的焦點。

以往，層合板的振動問題求解主要是基于理論解[1-2]。然而，由于數(shù)學(xué)問題的復(fù)雜性以及對實際問題進行建模的高要求，這些常規(guī)方法難以同時滿足計算效率和計算精度要求。隨著研究的不斷深入，人們提出了一系列的數(shù)值求解方法，如Galerkin法[3]、里茲法[4-6]、改進傅里葉級數(shù)法[7]、無網(wǎng)格法[8-15]等。Shi等用Galerkin法對四邊固定的任意層合板進行了自由振動分析。Ye等用Raleigh-Ritz法得到了具有一般邊界約束和內(nèi)線支承的復(fù)合材料中厚層板自由振動問題的修正傅立葉解。Shi等用改進的傅里葉級數(shù)方法(IFSM)研究了彈性地基上均勻支承和多點支承的中厚復(fù)合材料層合矩形板的自由和強迫振動特性。Xiang等提出了一種基于樣條徑向基函數(shù)和高階剪切變形理論的無網(wǎng)格法來分析固支層合板的自由振動問題。張科使用基于移動最小二乘(MLS-L)近似函數(shù)構(gòu)造的三維無網(wǎng)格法，對正交層合板的振動問題進行了分析。陳莘莘等提出了無網(wǎng)格自然鄰接點Pertov-Galerkin法分析了復(fù)合材料層合板的自由振動。使用插值型重構(gòu)核粒子無網(wǎng)格法，李情等分析了復(fù)合材料層合板的自由振動特性。通過對上述文獻分析可知，目前對于分析層合板振動特性的數(shù)值方法種類豐富，不同方法之間的區(qū)別主要是基/試函數(shù)和求解方法的不同，但是目前絕大多數(shù)的數(shù)值方法缺陷是在計算層合板振動問題建模過程復(fù)雜，并且求解過程一般采用傳統(tǒng)的數(shù)值積分方法，影響計算效率。

譜切比雪夫法(ST)采用具有指數(shù)收斂特性、優(yōu)越的數(shù)值穩(wěn)定性和高計算精度[16-17]的切比雪夫多項式作為基函數(shù)，來分析結(jié)構(gòu)的動力學(xué)問題。在文獻[18]中，Yagci等首次系統(tǒng)闡述了譜切比雪夫法，其核心思想在于：一是采用準(zhǔn)確的導(dǎo)數(shù)和內(nèi)積矩陣[19]來求解方程中導(dǎo)數(shù)和積分運算，以提高算法的計算效率；二是采用Gauss-Lobatto節(jié)點[20]對函數(shù)進行離散化處理以保證計算結(jié)果的準(zhǔn)確性?；谝陨蟽?yōu)點，譜切比雪夫法已成功應(yīng)用于各種結(jié)構(gòu)[21-26]的動力學(xué)特性研究中，例如Bediz使用二維譜切比雪夫(2D-ST)法預(yù)測了任意幾何形狀厚板在不同邊界條件下的動力特性，然而對于復(fù)合材料層合板的自由振動問題暫無涉及。

本文將譜切比雪夫法應(yīng)用于復(fù)合材料層合板的自由振動特性研究?；谝浑A剪切變形理論，引入邊界彈簧技術(shù)實現(xiàn)不同邊界條件的模擬，得到層合板的能量方程，應(yīng)用二維譜切比雪夫法求解能量方程獲得結(jié)構(gòu)的振動分析模型。通過數(shù)值分析算例，與其他方法的計算結(jié)果進行比較，證明了此方法的有效性，還在此基礎(chǔ)上研究了彈性模量比和鋪設(shè)角對復(fù)合材料層合板自由振動特性的影響。

1 理論公式

1.1 模型描述

如圖1所示為復(fù)合材料層合板的幾何形狀和坐標(biāo)系。直角坐標(biāo)系o-xyz位于結(jié)構(gòu)的幾何中面上，層合板沿x,y和z方向的長度、寬度和厚度分別是a、b和h；α是層纖維方向和x軸之間的鋪設(shè)角度；U、V、W分別表示層合板中任意一點在x、y、z方向上的位移。在本文中，采用邊界彈簧技術(shù)實現(xiàn)不同類型的邊界條件，引入了三組線性約束彈簧(ku、kv、kw)和兩組轉(zhuǎn)角約束彈簧(kφ、kθ)，kt(t=u,v,w,φ,θ)的值可以從零變化到無窮大，它們的不同組合可以表示各種邊界約束條件。例如，kt= 0(t=u,v,w,φ,θ)表示自由邊界條件(F)，而kt= ∞(t=u,v,w,φ,θ)表示固支邊界條件(C)。

圖1 復(fù)合材料層合板的幾何模型與坐標(biāo)系Fig.1 Geometric model and coordinate system of composite laminated plates

1.2 能量方程

應(yīng)用一階剪切變形理論，復(fù)合材料層合板中任意一點的位移場的表達式為

(1)

式中：廣義位移向量u= (uvwφθ)T，u,v和w表示在層合板幾何中面上沿著x,y和z方向的位移分量；φ和θ分別表示分別為x-z平面和y-z平面的轉(zhuǎn)角位移。

由式(1)得應(yīng)變-位移線性關(guān)系(幾何方程)為

(2)

(3)

根據(jù)廣義胡克定律，層合板第k層的應(yīng)力-應(yīng)變關(guān)系(物理方程)為

(4)

力和力矩的合力是通過積分板厚度上的應(yīng)力獲得的，通過引入剪切校正因子κ，最終可以得出廣義力與應(yīng)變之間的關(guān)系。其關(guān)系式表示如下

(5)

式中：N={NxNyNxy}T,M={MxMyMxy}T,QS={QxQy}T分別為合力，合力矩和等效剪力向量；A為拉伸剛度矩陣，B為耦合剛度矩陣，D為彎曲剛度矩陣，AS為剪切剛度矩陣，它們的表達式如下所示

(6)

因此，存儲在板中的應(yīng)變能US為

(7)

式中，R為5×5線性對稱微分算子矩陣。其表達式為

(8)

動能T的計算表達式為

(9)

利用方程(1)定義的位移場，沿z方向求積分，動能項可改寫為

(10)

G為5×5對稱質(zhì)量矩陣，其表達式為

(11)

式中：(I0,I1,I2)為慣性矩；ρk為第k層材料密度。Usprings是具有彈性均勻邊界條件的板的四邊的勢能，它是用五種彈簧均勻分布在四個邊來模擬。

(12)

式中,k=diag(kukvkwkφkθ)。

因此層合板的拉格朗日能量方程可寫為

L=T-US-Usprings

(13)

1.3 譜切比雪夫近似原理

切比雪夫多項式為遞歸和正交多項式，在區(qū)間x∈[-1,1]上，其表達式為

Cl(x)=cos(larc cos(x)),l=0,1,2,…,

(14)

應(yīng)用區(qū)間變換，對于任意區(qū)間[l1,l2]的函數(shù)y(x)，可以使用有限項的切比雪夫多項式展開表示為

(15)

在離散化系統(tǒng)(即確定采樣點)時，采用Gauss-Lobatto采樣方法[20]

l=0,1,…,M-1

(16)

擴展系數(shù)al與采樣點yl=y(xl) 之間存在一對一的映射關(guān)系，采樣函數(shù)y={y0,y1, …,yM-1}的向量與(截斷的)切比雪夫展開系數(shù)向量a={a0,a1, …,aM-1}之間的關(guān)系可以寫成

y=ΓBa

a=ΓFy

(17)

(18)

b=a

(19)

(20)

其中s= (l2-l1)/2為區(qū)間縮放系數(shù)。

應(yīng)用切比雪夫多項式的積分性質(zhì)，y(x)的積分表達式為

(21)

為了得到內(nèi)積矩陣，令h(x) =f(x)g(x)，假設(shè)fd=diag[f1f2…fM]，跟據(jù)式(17)和(21)有

(22)

令vTΓF= [V1…VM]，故式(22)可改寫為

(23)

對于二維問題，函數(shù)f(x,y)將產(chǎn)生第二階張量，其組成元素表達式為

flm=f(xl,ym)

(24)

式中，l= 1, …,Mx和m= 1, …,My是張量的指數(shù)。fkl為組成Mx×My二維采樣函數(shù)矩陣F的元素值。使用以下映射算法獲得采樣二維函數(shù)的矢量形式，為

fc=flm,c=(l-1)My+m

(25)

應(yīng)用二維譜切比雪夫法，在坐標(biāo)(x,y)中的每個位移點都用雙切比雪夫多項式展開表示為

(26)

式中:alm是雙展開系數(shù);Mx和My是在每個方向上使用的切比雪夫多項式數(shù)。將式(25)張量到矢量的映射應(yīng)用于展開系數(shù)，式(26)可改寫為

(27)

二維問題的展開系數(shù)與采樣函數(shù)值的關(guān)系為

f=ΨBa

a=ΨFf

(28)

(29)

對于微分運算，需要獲得二維問題中關(guān)于x和y方向的導(dǎo)數(shù)矩陣Dx和Dy。按照文獻[25]的方法來擴展導(dǎo)數(shù)矩陣，就可以分別得到相對于x和y的擴展導(dǎo)數(shù)矩陣Dx和Dy。使用向前和向后變換矩陣以及擴展的導(dǎo)數(shù)矩陣，有如下定義：

fx=ΨBbx=ΨBDxa

fy=ΨBby=ΨBDya

(30)

對于面積積分，使用內(nèi)積矩陣方法

?Af(x,y)g(x,y)dxdy=fTVg

(31)

V為內(nèi)積矩陣，其表達式為

V=x?y

(32)

1.4 振動求解

應(yīng)用二維譜切比雪夫法求解式(13)可得到

(33)

式中，M,KU和KS分別代表質(zhì)量矩陣。應(yīng)變能剛度矩陣和邊界彈簧剛度矩陣，其具體表達式如下

(34)

t=u,v,w,φ,θ

Vx=x?INy×Ny,Vy=INx×Nx?y

(35)

I為單位矩陣。

對于式(33)，采用弱形式的求解方法得到對未知變量求得偏導(dǎo)數(shù)方程

(36)

將式(36)以矩陣的形式表達，則結(jié)構(gòu)的振動特性方程可以寫為以下形式

(KU+KS-Mω2)A=0

(37)

式中:A代表未知系數(shù)的矢量;ω表示結(jié)構(gòu)的固有頻率。通過直接求解式(37)，可以得到板的頻率特性和相應(yīng)的特征向量。

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

基于前面得到的自由振動特征求解方程，本部分將對復(fù)合材料層合板的自由振動特性進行分析。為簡化研究，如無特別說明，在接下來的計算分析中，層合板的幾何參數(shù)為：a= 1 m,b= 1 m,h= 0.1 m；材料參數(shù)為：E2= 10 GPa、E1= 25E2、G12= 0.5E2、G13=G12、G23=0.2E2、μ12= 0.25、μ12=μ21、ρ= 1 500 kg/m3，κ= 5/6。如前所述，對于不同的邊界可以設(shè)置與之對應(yīng)的彈簧約束剛度實現(xiàn)，符號C、S及F分別用來表示固支、簡支、自由邊界。特別地，本文所涉及的邊界組合按逆時針方向設(shè)定。此外，本文方法具有簡便性和通用性，當(dāng)結(jié)構(gòu)的幾何參數(shù)和邊界條件發(fā)生變化時，只需要修改相關(guān)參數(shù)即可，并不用重新推導(dǎo)，從而可以大大簡化計算過程。此外，本文中采用式(38)對層合板的固有頻率進行無量綱化。

(38)

2.1 數(shù)值驗證

在本節(jié)中，通過與文獻中的結(jié)果進行比較，對二維譜切比雪夫法進行了數(shù)值驗證。首先，表1為改變x和y方向上所需的多項式數(shù)Mx×My，四端固支的方形層合板的前五階無量綱頻率。層合材料的纖維鋪設(shè)角為[30°/-30°]，板厚h=0.05 m。從表1可以看出，隨著多項式數(shù)的遞增，各階頻率依次從低階到高階收斂至穩(wěn)定值, 當(dāng)Mx×My=16×16時，所有無量綱頻率均已達到收斂值。與Jacobi-Ritz法進行對比，兩種方法均是基于一階剪切變形理論，所以兩者前五階無量綱頻率計算結(jié)果高度吻合，相對誤差在0.001%以內(nèi)。此外，本文的譜切比雪夫法達到收斂值的矩陣維數(shù)為1 280×1 280 (5MxMy× 5MxMy)，而文獻[6]的矩陣維數(shù)為5 120×5 120。本文只需要少量的多項式數(shù)即可獲得足夠準(zhǔn)確的結(jié)果，充分說明了本文所提出的譜切比雪夫法具有收斂速度快，同時計算準(zhǔn)確性好的特點。

表1 采用不同多項式數(shù) Mx×My的復(fù)合材料層合板的前五階無量綱頻率對比Tab.1 Comparison of the first five-order dimensionless frequency of composite laminated plates with different Mx×My

表2所示為在固支邊界條件下，[0/α]和[α/-α]兩種鋪設(shè)方式下的復(fù)合材料層合板的前三階無量綱頻率，密度為ρ=1 000 kg/m3。將計算結(jié)果分別和Galerkin法和局部移動Kriging無網(wǎng)格法的計算結(jié)果進行了比較,從表2可以看出，三種方法的計算結(jié)果可以很好地匹配，微小的誤差來源于主要是因為不同數(shù)值方法，不同方法之間的計算原理在本質(zhì)上存在一定差異。

表2 [0/α]和[α/-α]兩種鋪設(shè)方式下復(fù)合材料層合板的前三階無量綱頻率對比Tab.2 Comparison of the first three-order dimensionless frequencies of composite laminated plates with the two layer types [0/α] and [α/-α]

最后，由于本文方法是為了求解層合板在任意邊界條件下的自由振動問題，因此有必要驗證其在各種邊界條件下的有效性和計算精度，相關(guān)對比數(shù)據(jù)結(jié)果如表3所示，參考值為修正傅里葉解計算結(jié)果。材料密度ρ=1 000 kg/m3，纖維鋪設(shè)角為[45°/-45°]。在該算例中，考慮了四種邊界條件，分別是C-C-C-C、S-S-S-S、F-C-F-C、F-S-F-S，和兩種寬厚比，即h/a=0.05和0.2。兩種方法同時基于一階剪切變形理論，只是基函數(shù)和求解形式不同，所以計算結(jié)果非常接近。結(jié)果表明，本文方法適用于任意邊界條件。

表3 不同邊界條件下[45°/-45°]復(fù)合材料層合板的前三階無量綱頻率對比Tab.3 Comparison of the first three-order dimensionless frequencies of [45°/-45°] composite laminated plates with different boundary conditions

2.2 參數(shù)化研究

根據(jù)引言調(diào)研的文獻可知，在工程實踐中，板結(jié)構(gòu)的材料屬性往往是多種多樣的，而纖維的敷設(shè)角度和纖維自身的彈性模量取值對于層合板的振動特性而言至關(guān)重要，這也直接反應(yīng)在復(fù)合材料的纖維類型選擇和加工工藝上。因此，詳細(xì)研究纖維的敷設(shè)角度和纖維自身的彈性模量取值對復(fù)合材料層合板的自由振動特性的影響對指導(dǎo)實際工程設(shè)計和應(yīng)用有著重要意義。本節(jié)以四端固支，鋪設(shè)方式為[α/-α]的層合板為例，研究了彈性模量比E1/E2和鋪設(shè)角α對復(fù)合材料層合板無量綱頻率的影響。層合板的其余結(jié)構(gòu)參數(shù)參考默認(rèn)。E1/E2的變化區(qū)間為[1,60]，α從0°遞增至180°，提取層合板的1階、10階和100階無量綱頻率，最終計算結(jié)果如圖2所示。需要注意的是，當(dāng)計算100階固有頻率時，所采用的多項式數(shù)為Mx×My=40×40。從圖2可以看出，隨著E1/E2的增大，1階、10階和100階的無量綱頻率也隨之遞增，E1/E2取最大值60時，頻率也達到了最大值；對于鋪設(shè)角α，1階、10階和100階的無量綱頻率變化規(guī)律關(guān)于α=90°對稱，并且模態(tài)階次不同，頻率變化規(guī)律也會略有不同，1階頻率在0°、90°和180°附近達到最大值，10階和100階頻率則在區(qū)間[30°,60°]和[120°,150°]內(nèi)取最大值。

1階

10階

100階圖2 彈性模量比E1/ E2和鋪設(shè)角α對層合板無量綱頻率的影響Fig.2 Effects of elastic modulus ratio E1/ E2 and laying angle α on non-dimensional natural frequencies of laminated plates

3 結(jié) 論

本文在一階剪切變形理論的基礎(chǔ)上，應(yīng)用二維譜切比雪夫法建立了任意邊界條件下復(fù)合材料層合板的自由振動分析模型。通過數(shù)值對比分析和參數(shù)化研究，得出以下結(jié)論：

(1) 所采用的譜切比雪夫法具有收斂速度快，計算準(zhǔn)確性好的特點，可以有效分析任意邊界條件下復(fù)合材料層合板的自由振動問題。

(2) 彈性模量比E1/E2對復(fù)合材料層合板的自由振動頻率有重要影響。隨著彈性模量比E1/E2的增加，結(jié)構(gòu)的無量綱固有頻率也會顯著增加。

(3) 對于鋪設(shè)角α，復(fù)合材料層合板的無量綱頻率關(guān)于α=90°對稱分布，并且與模態(tài)階次有一定的關(guān)聯(lián)性。