杜 圓, 李海超, 龐福振, 繆旭弘,2
(1.哈爾濱工程大學(xué) 船舶工程學(xué)院,哈爾濱 150001; 2.海軍研究院,北京 100161)
矩形板結(jié)構(gòu)在船舶、海洋平臺(tái)、潛艇及波浪能電站等結(jié)構(gòu)中有著大量的應(yīng)用,因而設(shè)計(jì)者需對(duì)其振動(dòng)特性進(jìn)行詳細(xì)了解。目前常用的方法有能量法[1]、微分求積法[2]、辛本征展開法[3]。但使用能量法的前提是選取與邊界條件相對(duì)應(yīng)的撓度函數(shù),微分求積法中加權(quán)系數(shù)的選擇對(duì)結(jié)果影響很大,辛本展開法所建立的泛函缺乏與邊界條件相關(guān)的表征項(xiàng)。Shen等[4]通過使用坐標(biāo)函數(shù)的線性組合來表征彎曲位移函數(shù)得到了板結(jié)構(gòu)的彎曲振動(dòng)近似解;Cho等[5]基于假定模態(tài)法,通過對(duì)板的拉格朗日方程求解特征值求解板的自由振動(dòng)問題;Reddy[6]基于Reissner-Midlin剪切變形理論與非線性應(yīng)變位移關(guān)系對(duì)板的彎曲振動(dòng)進(jìn)行求解;Aksu等[7]運(yùn)用有限差分方法對(duì)矩形板固有頻率進(jìn)行了計(jì)算;Gorman[8]提出用系統(tǒng)疊加法來求解各種邊界條件下板的振動(dòng)問題。上述研究求解過程復(fù)雜,且僅能對(duì)經(jīng)典邊界條件求解,無法對(duì)實(shí)際工程應(yīng)用中常見的彈性邊界條件進(jìn)行求解。
本文將薄板振動(dòng)的位移容許函數(shù)表示成二維傅里葉余弦級(jí)數(shù)和輔助傅里葉級(jí)數(shù)的線性組合,輔助傅里葉級(jí)數(shù)的引入使位移容許函數(shù)三階連續(xù)且可偏導(dǎo);從而滿足薄板結(jié)構(gòu)振動(dòng)控制微分方程,解決了傳統(tǒng)傅里葉級(jí)數(shù)法在邊界處不連續(xù)的問題。基于改進(jìn)后的位移容許函數(shù)列出矩形薄板拉格朗日能量泛函方程,然后通過Hamilton原理求解得到位移容許函數(shù)中的系數(shù)。通過改變橫向位移約束彈簧剛度值k和旋轉(zhuǎn)約束彈簧剛度值K模擬任意邊界條件。
在與文獻(xiàn)及數(shù)值仿真結(jié)果進(jìn)行對(duì)比驗(yàn)證本文方法的有效性后,對(duì)彈性組合邊界條件下矩形薄板的自由振動(dòng)特性進(jìn)行求解。并通過大量計(jì)算探討邊界條件與矩形薄板自由振動(dòng)特性的關(guān)系,旨在為后續(xù)研究及工程應(yīng)用提供參考。
本文所研究的任意邊界條件下矩形薄板示意圖如圖1所示。矩形薄板的邊界條件通過四條邊上的橫向位移約束彈簧與旋轉(zhuǎn)約束彈簧來模擬。
圖1 任意邊界條件下的矩形薄板示意圖
橫向位移約束彈簧與旋轉(zhuǎn)約束彈簧在四條邊上均勻分布,其中kx0、ky0、kxa、kyb分別代表x=0,y=0,x=a,y=b邊界上橫向位移約束彈簧的剛度,Kx0、Ky0、Kxb、Kya分別代表x=0,y=0,x=a,y=b邊界上旋轉(zhuǎn)約束彈簧的剛度。當(dāng)kx0、ky0、kxa、kyb、Kx0、Ky0、Kxa、Kyb均趨于無窮時(shí),表征矩形薄板四邊剛性固定;當(dāng)kx0、ky0、kxa、kyb、Kx0、Ky0、Kxa、Kyb均設(shè)置為零時(shí),表征矩形薄板四邊自由;當(dāng)kx0、ky0、kxa、kyb剛度設(shè)置為無窮,Kx0、Ky0、Kxa、Kyb設(shè)置為0時(shí),表征矩形薄板四邊簡(jiǎn)支。
位移容許函數(shù)的選擇將直接影響后續(xù)計(jì)算精度[9],求解過程中需將式(1)代入結(jié)構(gòu)控制方程式(3)。但傳統(tǒng)傅里葉級(jí)數(shù)的導(dǎo)數(shù)在端點(diǎn)處仍會(huì)產(chǎn)生不連續(xù)或跳躍,如圖2所示。通過將位移容許函數(shù)表征為式(1)的形式,輔助函數(shù)的引入使其滿足求解域內(nèi)三階求導(dǎo)連續(xù)且四階導(dǎo)數(shù)存在,有效克服邊界不連續(xù)現(xiàn)象,如圖3所示。
式中:λm=mπ/a;λn=nπ/b;Amn、Bmn、Cmn和Dmn為矩形薄板位移容許函數(shù)中的Fourier系數(shù),簡(jiǎn)諧時(shí)間因子eiωt表示薄板位移隨時(shí)間的變化。式(1)中第四項(xiàng)傅里葉正弦級(jí)數(shù)之和用來處理邊界不連續(xù),前三項(xiàng)中的m與n取值越大,位移容許函數(shù)越精確;在計(jì)算過程中展開級(jí)數(shù)項(xiàng)m與n的取值M、N稱為截?cái)嘀怠?/p>
圖2 傳統(tǒng)傅里葉級(jí)數(shù)在端點(diǎn)位置可能產(chǎn)生不連續(xù)問題示意圖
圖3 改進(jìn)傅里葉級(jí)數(shù)如何解決端點(diǎn)位置可能產(chǎn)生不連續(xù)問題示意圖
薄板僅承受法向載荷q(x,y,t)時(shí),薄板運(yùn)動(dòng)微分方程如式(2)所示
(2)
本文僅研究矩形薄板的自由振動(dòng),薄板上并不受力,結(jié)構(gòu)控制方程表示為
D▽4w(x,y)-ρhω2w(x,y)=0
(3)
式中:D為薄板的彎曲剛度,表達(dá)式為D=Eh3/(12(1-μ2));E為薄板的楊氏模量;μ為薄板的泊松比;ρ為薄板的密度;h為薄板的厚度;w(x,y)為“1.2”中的位移容許函數(shù)。
矩形薄板彎曲產(chǎn)生的應(yīng)變能如式(4)所示
(4)
矩形薄板邊界彈簧所儲(chǔ)存的彈性勢(shì)能如式(5)所示
(5)
矩形薄板系統(tǒng)的動(dòng)能T表示為
(6)
綜上,將式(4)到式(6)代入可得矩形薄板拉格朗日函數(shù)為
L=Vp+Vs-T
(7)
由Hamilton原理可知
(8)
將式(7)代入式(8)求解可得矩陣方程如下
([K]-ω2[M]){B}=0
(9)
式中:[K]=[Kp]+[Ks];[Kp]為矩形薄板應(yīng)變勢(shì)能剛度矩陣;[Ks]為矩形薄板邊界彈簧勢(shì)能剛度矩陣;M為矩形薄板質(zhì)量矩陣;B為未知的傅里葉級(jí)數(shù)矩陣;式(9)中各項(xiàng)可參照附錄,對(duì)式(9)求解可得到矩形薄板的固有頻率及振型。
為驗(yàn)證本文方法有效性,以下計(jì)算中矩形薄板材料參數(shù)參照文獻(xiàn)[10]設(shè)置如下:密度ρ=7 850 kg/m3,泊松比μ=0.3,楊氏模量E=2.1×1011Pa;為了便于與已有文獻(xiàn)[10]及有限元計(jì)算結(jié)果進(jìn)行對(duì)比,對(duì)計(jì)算結(jié)果進(jìn)行無量綱化處理,無量綱化公式為:Ω=ωa2(ρh/D)1/2。有限元計(jì)算通過ABAQUS實(shí)現(xiàn),材料參數(shù)與文獻(xiàn)[10]一致,單元類型為S4R,單元數(shù)為105。為便于表示,下文中邊界條件:C為剛固,S為簡(jiǎn)支,F(xiàn)為自由,E為彈性。
由上述分析可知,矩形薄板邊界條件通過四條邊上均布的橫向位移約束彈簧與旋轉(zhuǎn)約束彈簧來模擬;現(xiàn)討論彈簧剛度值對(duì)矩形薄板固有頻率的影響,矩形薄板邊長(zhǎng)a=1 m,寬b=1 m,厚度h=0.002 m。簡(jiǎn)支邊界條件下逐漸增大旋轉(zhuǎn)約束彈簧剛度,矩形薄板無量綱頻率參數(shù)變化,如圖4所示。旋轉(zhuǎn)約束彈簧剛度大于1010N·m/rad后,矩形薄板無量綱頻率參數(shù)基本不變;旋轉(zhuǎn)約束彈簧剛度設(shè)為1010N·m/rad時(shí),矩形薄板無量綱頻率參數(shù)隨橫向位移約束彈簧剛度變化,如圖4所示。
圖4 矩形薄板無量綱頻率參數(shù)隨旋轉(zhuǎn)約束彈簧剛度值變化關(guān)系
圖5 矩形薄板無量綱頻率參數(shù)隨橫向位移彈簧剛度值變化關(guān)系
由圖4及圖5可知,當(dāng)旋轉(zhuǎn)約束彈簧剛度與橫向位移約束彈簧剛度大于1010時(shí),矩形薄板無量綱頻率參數(shù)已收斂,可作為剛性邊界條件;約束彈簧剛度均設(shè)置為0時(shí),可視為自由邊界條件;橫向位移約束彈簧與旋轉(zhuǎn)約束彈簧剛度值分別取102~107N/m與102~105N·m/rad的某一值時(shí),該邊界條件介于自由邊界條件與剛性固定邊界條件之間,可視為某種彈性邊界條件。綜上所述,參照文獻(xiàn)[11]本文任意邊界彈簧剛度取值如表1所示。
表1 任意邊界彈簧剛度值
截?cái)嘀礛、N取值不同,會(huì)對(duì)計(jì)算精度與效率產(chǎn)生較大影響;因而有必要對(duì)截?cái)嘀蹬c本文計(jì)算方法收斂性之間關(guān)系進(jìn)行探究。以S-S-S-S(四邊簡(jiǎn)支)邊界條件下的矩形薄板為研究對(duì)象,矩形薄板邊長(zhǎng)a=1 m,寬b=1 m,厚度h=0.002 m。簡(jiǎn)支邊界條件中橫向位移約束彈簧剛度k為1010N/m,旋轉(zhuǎn)約束彈簧剛度K設(shè)置為0。
表2為截?cái)嘀礛=N取值從5變?yōu)?4時(shí),該矩形薄板的前5階無量綱頻率參數(shù),此外有限元分析結(jié)果與參考文獻(xiàn)[10]結(jié)果也列在表2中作為參考。由表2可知,截?cái)鄶?shù)M=N=13后,算例中矩形薄板固有頻率已基本不變,本文計(jì)算方法已收斂,在后面的計(jì)算中將選取截?cái)鄶?shù)M=N=13。
表2不同的截?cái)嘀礛、N下S-S-S-S矩形薄板結(jié)構(gòu)無量綱頻率參數(shù)Ω
Tab.2 Dimensionless frequency parameters of S-S-S-S boundary condition with different truncated numbers
M=N模態(tài)階次12345521.67251.30651.31080.919100.676621.67151.30551.30880.918100.674721.67151.30351.30580.918100.667821.67151.30251.30380.917100.666921.67051.30251.30380.917100.6641021.67051.30251.30280.917100.6631121.67051.30151.30180.917100.6621221.67051.30151.30180.917100.6611321.67051.30151.30180.917100.6611421.67051.30151.30180.917100.661FEM21.65051.28051.29080.890100.650文獻(xiàn)[10]21.60051.28751.28780.886100.680
對(duì)應(yīng)表2中任意給定階次的模態(tài)頻率,式(9)均可求解得到相對(duì)應(yīng)的特征向量;特征向量對(duì)應(yīng)該矩形板在給定階次模態(tài)頻率所對(duì)應(yīng)振型的傅里葉展開系數(shù),將系數(shù)代入式(1)所表征的位移容許函數(shù),即可得到與給定階次模態(tài)頻率相對(duì)應(yīng)的模態(tài)陣型。圖6給出了S-S-S-S邊界條件下采用本文計(jì)算方法得到的矩形薄板部分模態(tài)振型圖。使用有限元分析軟件(ABAQUS)計(jì)算,單元類型為S4R,單元數(shù)為105時(shí)得到的部分模態(tài)振型圖如圖7所示。
(a) 第1階模態(tài)振型
(b) 第4階模態(tài)振型
(a) 第1階模態(tài)振型
(b) 第4階模態(tài)振型
表3~表5列出S-S-F-F、S-S-S-S、C-F-F-F三種典型邊界條件下,厚度h=0.002 m,不同長(zhǎng)寬比矩形薄板使用本文方法與有限元法計(jì)算得到的無量綱頻率參數(shù)及文獻(xiàn)解。
表3 S-S-F-F邊界條件下不同長(zhǎng)寬比矩形薄板無量綱頻率參數(shù)Ω
表4 S-S-S-S邊界條件下不同長(zhǎng)寬比矩形無量綱頻率參數(shù)Ω
表5 C-F-F-F邊界條件下不同長(zhǎng)寬比矩形無量綱頻率參數(shù)Ω
表3~表5中采用本文方法得到的解與有限元及文獻(xiàn)解吻合良好。此外,圖6~圖7中采用本文方法與有限元得到的振型圖趨勢(shì)相仿;以上對(duì)比說明,本文方法準(zhǔn)確可靠。此外,基于本文方法MATLAB求解時(shí)間僅需3 s(PC、3.4 GHz),而有限元求解超過60 s;本方法對(duì)不同尺寸、邊界條件的矩形薄板求解,僅需修改參數(shù)無需重新建模。
基于上述研究?jī)?nèi)容,對(duì)任意邊界條件下不同長(zhǎng)寬比,厚h=0.002 m的矩形薄板自由振動(dòng)特性進(jìn)行分析,旨在為后續(xù)研究及工程應(yīng)用提供參考。從表6及表7可知,邊界條件對(duì)矩形薄板無量綱頻率參數(shù)有顯著影響。固支邊界條件與彈性邊界條件組合中,隨著固支邊條界的增多,矩形薄板無量綱頻率參數(shù)呈增大趨勢(shì)。簡(jiǎn)支及自由邊界條件與彈性邊界條件組合中,隨著彈性邊條界的增多,矩形薄板無量綱頻率參數(shù)呈增大趨勢(shì)。
表6 不同邊界矩形薄板自由振動(dòng)無量綱頻率參數(shù)Ω
Tab.6 Dimensionless frequency parameters of rectangular thin plates in different boundary conditions
長(zhǎng)寬比r模態(tài)階次邊界條件C-C-C-EC-C-E-EC-E-E-EE-E-E-ES-S-S-E133.75432.18930.45729.16922.182259.26352.22550.41749.19248.237369.74767.68555.33349.19648.6041491.67877.00070.60564.53271.725596.48983.14374.25970.90382.6646127.065108.36290.44774.59995.727156.41453.15952.53651.98340.135287.37883.08178.90175.67664.5893118.712103.699102.989102.397100.2241.54141.000127.412119.885110.025109.3605144.928135.956122.755119.418121.1156193.532159.885156.656146.593161.383191.25485.86685.53385.24365.6622117.624110.700108.748107.02688.7293166.273157.983150.609145.346131.27224204.100178.161177.750177.380173.2815226.553197.675195.366193.435192.9486237.579228.259208.456194.921194.4181137.334129.299129.079128.89998.6302160.418150.510149.353148.317120.7813204.113192.221188.439185.107161.4212.54270.414256.904246.754238.450222.2555314.536274.658274.361274.119267.2846335.077291.831290.393289.136286.3271194.204182.963182.807182.686138.9972215.099201.769201.019200.286160.5283254.715238.789236.491234.081199.64334315.980297.729291.730286.417258.3575400.076379.768366.959354.642337.5436449.806392.928392.731392.503382.2521252.189246.701246.619246.488186.7062274.489263.862263.312262.773207.8513281.078297.287295.624294.112245.7633.54351.198350.970347.072343.349302.6085374.110427.171418.886411.329379.5556502.399526.706510.274495.878477.109
通過將薄板振動(dòng)的位移函數(shù)表示成二維傅里葉余弦級(jí)數(shù)和輔助傅里葉級(jí)數(shù)的線性組合,列出矩形薄板的拉格朗日方程;然后通過Hamilton原理對(duì)拉格朗日方程求解得到矩形薄板自由振動(dòng)無量綱頻率參數(shù)。通過與文獻(xiàn)及有限元計(jì)算結(jié)果對(duì)比,驗(yàn)證了本文方法的有效性。通過本文研究,可得如下主要結(jié)論:
表7 不同邊界矩形薄板自由振動(dòng)無量綱頻率參數(shù)Ω
Tab.7 Dimensionless frequency parameters of rectangular thin plates in different boundary conditions
長(zhǎng)寬比r模態(tài)階次邊界條件S-S-E-ES-E-E-EF-F-F-EF-F-E-EF-E-E-E124.71426.7568.18519.54320.739247.56148.3588.44522.17331.896348.00948.74418.95135.13643.9581467.19965.54027.87042.92653.154573.37372.93328.06845.91055.835693.92480.17846.19958.42868.620148.65750.0717.74744.08945.157267.42171.59914.06546.93255.3413101.073101.68631.41458.62080.8311.54108.026109.20045.88986.97597.6975116.744118.13653.86396.713106.5496151.198148.48868.45899.770120.851182.88483.89421.08278.48679.494298.861102.69842.93981.36688.9553134.042139.90063.96992.678112.25024176.394176.84973.902117.583153.3335190.238191.81291.444162.818173.0056192.035193.679113.479172.061181.5251127.117127.89129.437122.742123.6822141.459144.62447.347125.594132.6833172.362178.58883.824136.761154.0172.54224.676231.849111.930159.639192.4285273.338273.708136.638200.164249.5906286.278287.633147.401262.685269.8921181.266181.8959.283176.935177.5812194.628197.23339.298179.648186.3753222.427228.08060.702191.268205.98734269.682278.07695.207212.410241.4575339.301347.347158.304251.670296.2256391.881392.188191.591306.672369.0811245.340245.83550.661240.835241.7252258.040260.22974.098243.526250.2443283.652288.590108.545254.648268.9653.54326.761334.866149.489275.451302.1475390.935401.180232.355310.261353.0296477.963487.258247.703362.905423.776
(1) 本文算法收斂性與邊界彈簧剛度及截?cái)鄶?shù)相關(guān)。旋轉(zhuǎn)約束彈簧剛度K>1010N·m/rad,橫向位移約束彈簧剛度k>1010N/m時(shí),矩形薄板無量綱頻率參數(shù)已收斂,可作為剛性邊界條件。橫向位移約束彈簧剛度k為102N/m~107N/m中某一值時(shí),可將其視為某種彈性邊界條件。
(2) 典型邊界條件下本文計(jì)算結(jié)果與文獻(xiàn)及有限元解吻合良好,方法準(zhǔn)確可靠。
(3) 本文方法與有限元相比具有更高的求解效率;對(duì)不同尺寸、邊界條件的矩形薄板求解,僅需修改參數(shù),無需重新建模。
(4) 固支邊界條件與彈性邊界條件組合中,隨著固支邊條界范圍增大,矩形薄板無量綱頻率參數(shù)呈增大趨勢(shì)。簡(jiǎn)支及自由邊界條件與彈性邊界條件組合中,隨著彈性邊條界范圍增大,矩形薄板無量綱頻率參數(shù)呈增大趨勢(shì)。
[cos(λm1x)cos(λn1y)]dxdy
[cos(λm1x)sin(λn1y)]dxdy
MBA=MABT,MCA=MACT
2(1-μ)λmλnλm1λn1sin(λmx)sin(λny)sin(λm1x)sin(λn1y)]dxdy
2(1-μ)λmλnλm1λn1sin(λmx)sin(λny)sin(λm1x)cos(λn1y)]dxdy