張 俊， 李天勻， 朱 翔， 郭文杰

(1.華中科技大學(xué) 船舶與海洋工程學(xué)院，武漢 430074； 2.船舶與海洋水動(dòng)力湖北省重點(diǎn)實(shí)驗(yàn)室，武漢 430074；3.高新船舶與深海開發(fā)裝備協(xié)同創(chuàng)新中心，上海 200240; 4.華東交通大學(xué) 鐵路環(huán)境振動(dòng)與噪聲教育部工程研究中心，南昌 330013)

板殼結(jié)構(gòu)廣泛用于工程領(lǐng)域中，薄板振動(dòng)分析的解析方法，大多只能針對(duì)簡單形狀，然而實(shí)際問題卻非常復(fù)雜，因此對(duì)復(fù)雜形狀板的振動(dòng)問題進(jìn)行研究具有重要意義。

對(duì)任意形狀問題求解主要采用數(shù)值方法，如有限元方法、微分求積法、邊界點(diǎn)法等。數(shù)值方法以其天然的優(yōu)勢較為容易地應(yīng)用于求解任意形狀結(jié)構(gòu)問題，但同時(shí)也存在一些劣勢。比如劃分網(wǎng)格依賴使用者的經(jīng)驗(yàn)，在計(jì)算精度分析上需要耗費(fèi)精力；模型參數(shù)有變動(dòng)時(shí)需重新建模，比較費(fèi)時(shí)費(fèi)力。

解析方法作為一種簡單有效的分析方法，已經(jīng)廣泛地應(yīng)用于簡單板殼問題的靜動(dòng)態(tài)分析。焦映厚等[1]采用了波動(dòng)法研究有限尺寸加肋L型板結(jié)構(gòu)的振動(dòng)特性。楊念等[2]應(yīng)用解析方法研究了含復(fù)雜結(jié)構(gòu)應(yīng)力的平板振動(dòng)問題。Rayleigh-Ritz法作為一種常見的解析方法，得到廣泛的應(yīng)用。李正良等[3]基于Rayleigh-Ritz法分析了正交加筋圓柱殼-球殼組合結(jié)構(gòu)的自振特性。邱永康等[4]應(yīng)用Rayleigh-Ritz法，對(duì)任意邊界條件下中心開口矩形板的自振特性進(jìn)行分析。李凱等[5]結(jié)合Rayleigh-Ritz法計(jì)算了開口矩形板的固有頻率和振型。解析方法具有較多的優(yōu)勢，能夠從理論上揭示結(jié)構(gòu)振動(dòng)特性，便于機(jī)理分析，可以用于檢驗(yàn)數(shù)值方法等。然而對(duì)于復(fù)雜形狀的平板，試函數(shù)很難得到，公式推導(dǎo)非常復(fù)雜等，傳統(tǒng)Rayleigh-Ritz法等解析方法難以解決復(fù)雜形狀板殼振動(dòng)問題。

解析方法要實(shí)現(xiàn)對(duì)復(fù)雜形狀問題的求解，對(duì)其帶來的曲邊的邊界條件問題的處理尤為關(guān)鍵。虛擬彈簧模擬邊界條件作為一種靈活的方法被廣泛采用[6-9]。Li[10]應(yīng)用解析法分析了矩形板的面內(nèi)振動(dòng)，將位移彈簧和轉(zhuǎn)角彈簧取相應(yīng)的剛度值來模擬邊界條件。但是以上研究多是處理直邊邊界，對(duì)彈簧模擬曲邊的邊界條件的研究還較少。

在Rayleigh-Ritz法的基礎(chǔ)上，本文提出了一種用于求解復(fù)雜形狀板振動(dòng)問題的半解析方法——改進(jìn)的Rayleigh-Ritz法。傳統(tǒng)的Rayleigh-Ritz法及目前提出的一些計(jì)算特定形狀平板振動(dòng)的方法適應(yīng)性不太好，且很多形狀無法找到合適的試函數(shù)。本文方法適用于計(jì)算復(fù)雜形狀平板振動(dòng)問題，在不重新建模的情況下，可以對(duì)各種復(fù)雜的邊界條件進(jìn)行計(jì)算，具有很高的效率和適應(yīng)性。方法主要改進(jìn)在：一是拓展試函數(shù)域，使模擬復(fù)雜形狀結(jié)構(gòu)的位移函數(shù)易得；二是結(jié)合彈簧邊界方法，用于處理復(fù)雜形狀結(jié)構(gòu)的邊界條件，尤其是分析曲形邊界。與傳統(tǒng)Rayleigh-Ritz法相比，本方法的應(yīng)用范圍更廣，精度也較高，為工程問題中復(fù)雜形狀板的振動(dòng)問題的求解提供了參考。

1 改進(jìn)的Rayleigh-Ritz法

復(fù)雜形狀板物理模型如圖1所示：板面為S，沿邊界分布著位移線彈簧和轉(zhuǎn)角線彈簧，兩種彈簧的剛度系數(shù)分別為k和K。

圖1 復(fù)雜形狀薄板的物理模型

1.1 方法的主要特點(diǎn)

1.1.1 函數(shù)域的拓展

對(duì)Rayleigh-Ritz法中的函數(shù)域進(jìn)行擴(kuò)展，模型如圖2所示。

圖2 復(fù)雜形狀板的計(jì)算模型

S為薄板的求解域，Sa為試函數(shù)的積分域，長度分別為Lx和Ly。由于矩形區(qū)域的試函數(shù)易于表達(dá)，當(dāng)求解域是不規(guī)則形狀時(shí)，可以假想將結(jié)構(gòu)域擴(kuò)展為一個(gè)略大于真實(shí)域的矩形域，矩形域要能夠包含結(jié)構(gòu)域，按照矩形域選取結(jié)構(gòu)的試函數(shù)。計(jì)算時(shí)對(duì)結(jié)構(gòu)域進(jìn)行積分，得到板的應(yīng)變能、動(dòng)能等，結(jié)果非常準(zhǔn)確。研究表明，函數(shù)域越接近結(jié)構(gòu)域，得到的解越精確，這也是本方法選擇函數(shù)域的一個(gè)原則。所以當(dāng)矩形域剛好與不規(guī)則結(jié)構(gòu)域外切時(shí)，得到的結(jié)果誤差最小。

在求解應(yīng)變能、動(dòng)能時(shí)需對(duì)板的結(jié)構(gòu)域積分，以半徑為R的圓板為例，由于其對(duì)稱性，計(jì)算積分時(shí)只計(jì)算其1/4的結(jié)構(gòu)，彎曲應(yīng)變能為

(1)

(2)

圓板的應(yīng)變能和動(dòng)能通過二重積分得到。通常情況下，不規(guī)則形狀的邊界可以用多項(xiàng)式擬合，可以通過積分計(jì)算其能量。

值得注意的是，有些復(fù)雜的函數(shù)很難直接積分，為了解決這個(gè)難題，同時(shí)也為了提高計(jì)算效率，本文引入離散方法計(jì)算積分，具體離散的方式如圖3所示。

圖3 沿縱坐標(biāo)離散示意圖

沿y方向?qū)⒎e分域分為Q等份，每份寬為Δy，縱坐標(biāo)為yj，根據(jù)曲線的表達(dá)式計(jì)算出yj處對(duì)應(yīng)積分域的長度Lxj，表示為：

Lxj=xR-xL

(3)

再將Lxj平均分為Q份，長度為

Δxj=Lxj/Q

(4)

離散后小分段的橫縱坐標(biāo)分別為

xij=Lxj+(i-0.5)·Δxj

(5a)

yj=ymin+(i-0.5)·Δy

(5b)

則二重積分可以轉(zhuǎn)化為對(duì)微面積的求和，應(yīng)變能與動(dòng)能都可以表示為如下形式

(6)

式中:C為系數(shù);f(xij,yi)為計(jì)算應(yīng)變能與動(dòng)能時(shí)的相關(guān)的函數(shù)。

這樣大大提高了計(jì)算效率，并解決了部分函數(shù)難以直接積分計(jì)算的難題。

1.1.2 虛擬彈簧處理曲形邊界條件

經(jīng)典邊界條件，如固支和簡支等往往不能準(zhǔn)確模擬實(shí)際結(jié)構(gòu)的邊界，考慮一般性，采用彈簧模型模擬邊界條件。在結(jié)構(gòu)域的邊界上采用位移線彈簧和轉(zhuǎn)角線彈簧。

假設(shè)位移線彈簧和轉(zhuǎn)角線彈簧剛度值分別為kij(N/m2)和Kij(N/rad)，復(fù)雜邊界條件可通過調(diào)節(jié)彈簧的剛度值模擬。經(jīng)典邊界條件的剛度取值如表1所示。模擬彈性邊界條件時(shí)則將剛度系數(shù)取為相應(yīng)的值。

表1 經(jīng)典邊界條件下的彈簧取值

對(duì)于本文所研究的復(fù)雜形狀板，其邊界為曲形邊界，則彈簧的彈性勢能表示為

(7)

圖4 邊界處偏導(dǎo)數(shù)計(jì)算示意圖

其中曲邊對(duì)法線方向的偏導(dǎo)數(shù)表示為

(8)

對(duì)邊界做離散處理，如圖5所示。

圖5 邊界處離散示意圖

每個(gè)微段的長度為

(9)

則邊界處彈性勢能可以表示為

(10)

1.2 分析過程

薄板的位移函數(shù)可以表示為

(11)

式中:Amn為未知的展開系數(shù);m和n為序號(hào);M和N為截?cái)囗?xiàng)數(shù)；fm(x)和gn(y)為沿x和y方向的正交多項(xiàng)式，選用改進(jìn)的傅里葉級(jí)數(shù)，eiωt為簡諧時(shí)間因子。

應(yīng)用本文方法時(shí)，選取合適的試函數(shù)至關(guān)重要。由于改進(jìn)的傅里葉級(jí)數(shù)可以滿足任意邊界條件，使板的試函數(shù)在整個(gè)求解域內(nèi)三階導(dǎo)數(shù)連續(xù)且四階導(dǎo)數(shù)各點(diǎn)均存在，有效地克服邊界處可能出現(xiàn)的不連續(xù)現(xiàn)象，選用改進(jìn)的傅里葉級(jí)數(shù)作為試函數(shù)[11]，表示為

fm(x)=sin(mπx) 0

fm(x)=cos[(m-5)πx]m≥5

(12a)

gn(x)=sin(nπx) 0

gn(x)=cos[(n-5)πx]n≥5

(12b)

式中:m=1，2，3，…,M；n=1，2，3，…,N。

為適用于計(jì)算復(fù)雜形狀平板在多種邊界條件下的振動(dòng)問題，選擇的試函數(shù)具有以下特點(diǎn)：是滿足相應(yīng)控制方程的正交函數(shù)；在邊界處要具有任意性，即不能天然的滿足某種幾何邊界條件。綜上，還可以選用切比雪夫級(jí)數(shù)[12]等符合上述條件的具有任意性的函數(shù)作為試函數(shù)。

薄板的彎曲應(yīng)變能可表示為

(13)

式中:w(x,y)為結(jié)構(gòu)振動(dòng)試函數(shù);μ為泊松比，D=Eh3/(12(1-μ2))為薄板的彎曲剛度，h為板厚，時(shí)間項(xiàng)eiωt可以解耦并省略。

忽略邊界彈簧的質(zhì)量，結(jié)構(gòu)的動(dòng)能為

(14)

式中:ρ為質(zhì)量密度。

系統(tǒng)的能量泛函可以表示為

∏=Vp+Vs-T

(15)

將式(7)、式(13)、式(14)代入式(15)，并對(duì)未知系數(shù)求極值：

(16)

式中:Amn是描述薄板彎曲振動(dòng)的未知系數(shù)。

于是結(jié)構(gòu)的振動(dòng)問題轉(zhuǎn)化成了求解特征值的問題，可以表示為

(K-ω2M)A=0

(17)

式中:K為剛度矩陣K=Ks+Kp，Ks為彈簧能量的剛度矩陣，Kp為整體結(jié)構(gòu)應(yīng)變能的剛度矩陣；M為結(jié)構(gòu)的質(zhì)量矩陣;A為未知的系數(shù)向量;ω為圓頻率。

2 數(shù)值分析

2.1 收斂性分析

試函數(shù)中的截?cái)囗?xiàng)數(shù)M、N的取值和彈簧剛度取無窮時(shí)的實(shí)際計(jì)算值對(duì)結(jié)果的精度有很大的影響，本次收斂性分析針對(duì)以上兩個(gè)量進(jìn)行。

選用四邊自由和固支的矩形板進(jìn)行收斂性分析，幾何參數(shù)如下：長a=2 m，寬b=1 m，厚度h=0.02 m。板的材料參數(shù)為：楊氏模量E=2.1×1011Pa，泊松比μ=0.3，密度ρ=7 850 kg/m3。

首先對(duì)截?cái)囗?xiàng)數(shù)進(jìn)行收斂性分析。表2給出剛度系數(shù)為0時(shí)，不同截?cái)囗?xiàng)數(shù)M=N=q時(shí)對(duì)應(yīng)的前10階固有頻率。

表2 矩形板固有頻率隨截?cái)囗?xiàng)數(shù)的變化

從表2可知，隨著截?cái)囗?xiàng)數(shù)M，N的增加，矩形板的自由振動(dòng)固有頻率趨于一致。當(dāng)M=N=12時(shí)，認(rèn)為算例已經(jīng)收斂，以下計(jì)算中取M=N=12。

對(duì)彈簧的剛度系數(shù)進(jìn)行收斂性分析。在處理固支等邊界時(shí)需將邊界彈簧的剛度系數(shù)取為無窮，實(shí)際計(jì)算時(shí)需要用一個(gè)較大的值來代替。表3為彈簧剛度系數(shù)取不同值時(shí)固有頻率的變化。算例中假設(shè)k和K均取同一數(shù)值，即k=pN/m2,K=pN/rad，D為板的彎曲剛度值。

表3 矩形板固有頻率隨彈簧剛度的變化

表3中數(shù)據(jù)表明：隨著剛度系數(shù)的逐漸增大，固有振動(dòng)頻率逐漸收斂。當(dāng)剛度系數(shù)K和k的數(shù)值大小均為108D時(shí)已經(jīng)收斂。故在彈簧邊界剛度為無窮大時(shí)，將彈簧剛度系數(shù)值取108D代入計(jì)算即可。

2.2 準(zhǔn)確性分析

本節(jié)驗(yàn)證方法計(jì)算直邊板、斜邊板、曲邊板及不規(guī)則形狀板自由振動(dòng)的準(zhǔn)確性，方法計(jì)算薄板問題效率較高，應(yīng)用MATLAB，得到較為精確的前10階固有頻率的時(shí)間僅需約40 s(PC，4 GHz，16 GB)。

2.2.1 矩形板的彎曲自由振動(dòng)分析

(18)

式中:fref為文獻(xiàn)結(jié)果或者有限元解。

表4 矩形板固有頻率

從表4可知，兩種結(jié)果對(duì)比誤差很小，說明本方法計(jì)算直邊板的自由振動(dòng)時(shí)的準(zhǔn)確性。

2.2.2 三角形板的彎曲自由振動(dòng)分析

為驗(yàn)證本文方法計(jì)算斜邊板的準(zhǔn)確性，對(duì)不同邊界條件下三角形板的彎曲自由振動(dòng)頻率進(jìn)行分析(物理模型見圖6)。

圖6 三角形板物理模型

2.2.3 橢圓形板的彎曲自由振動(dòng)分析

為驗(yàn)證本文方法求解曲邊板自由振動(dòng)的準(zhǔn)確性，下文將對(duì)在不同邊界條件下的橢圓形板自由振動(dòng)進(jìn)行分析(物理模型見圖8)。

圖7 簡支邊界條件下三角形形薄板高階振型與ansys軟件對(duì)比圖

表5 簡支三角形板的無量綱固有頻率

圖8 橢圓形薄板物理模型

橢圓板長軸a=4 m，短軸b=2 m板厚h=0.05 m。選取簡支和固支邊界條件對(duì)其進(jìn)行分析。

從表7可知，本方法計(jì)算結(jié)果與有限元(Ansys)結(jié)果吻合良好，最大誤差不超過2%。圖9給出了固支邊界下前四階模態(tài)振型圖，證明了本方法的準(zhǔn)確性。

2.2.4 曲邊梯形板的自由振動(dòng)分析

綜合上述幾種形狀，為了說明本方法求解復(fù)雜形狀平板自由振動(dòng)特性的準(zhǔn)確性，本節(jié)研究了帶有曲邊的梯形的自振特性，模型如圖10所示。其中“②”號(hào)和“④”號(hào)邊為曲邊，其參數(shù)如下：上底長a=4 m，下底長b=10 m，高c=8 m，兩斜邊為1/4橢圓曲線，其中長軸為8 m，短軸長度為3 m。用本方法計(jì)算該曲邊梯形板的固有頻率，并與Ansys軟件結(jié)果進(jìn)行對(duì)比。

表8給出前10階固有頻率結(jié)果對(duì)比。從表8可知，兩種方法結(jié)果吻合良好，誤差不超過1%。表明本文方法在計(jì)算曲邊直邊組合形狀時(shí)的準(zhǔn)確性。

表6 兩斜邊簡支底邊彈性邊界下三角形板固有頻率

(a)

(b)

(c)

(d)

表7 橢圓形薄板自由振動(dòng)固有頻率

表8 曲邊梯形板固有頻率

圖10 曲邊梯形板

3 結(jié) 論

本文基于改進(jìn)Rayleigh-Ritz法對(duì)復(fù)雜形狀薄板的自由振動(dòng)問題進(jìn)行分析。首先以矩形板為例進(jìn)行收斂性分析，再對(duì)矩形、三角形、橢圓形和曲邊梯形板的計(jì)算結(jié)果與有限元結(jié)果進(jìn)行對(duì)比，證明了本方法的準(zhǔn)確性。給出本文方法的低階和高階振型圖與有限元的振型圖進(jìn)行對(duì)比，直觀地看到振型吻合很好?？梢缘贸鼋Y(jié)論：改進(jìn)的Rayleigh-Ritz法在計(jì)算復(fù)雜形狀板的振動(dòng)問題上具有良好的收斂性、準(zhǔn)確性和適應(yīng)性。本文的方法從理論上可以拓展到對(duì)復(fù)雜形狀中厚板、厚板、功能梯度板、殼類動(dòng)態(tài)問題的求解。