廣義協(xié)調(diào)COONS曲面厚薄板通用單元振動(dòng)分析

陳　斌，張延慶

(北京工業(yè)大學(xué)建筑工程學(xué)院，北京　100124)

?

廣義協(xié)調(diào)COONS曲面厚薄板通用單元振動(dòng)分析

陳斌，張延慶

(北京工業(yè)大學(xué)建筑工程學(xué)院，北京100124)

摘要：厚薄板通用單元一直是板殼有限元研究中的重點(diǎn)和難點(diǎn)。根據(jù)單向厚板的基本方程推導(dǎo)得到內(nèi)含剪切變形因子的插值基函數(shù)，基于Mindlin-Reissner厚板理論并結(jié)合廣義協(xié)調(diào)思想，通過(guò)第2類(lèi)有限元COONS曲面構(gòu)造法及撓度和轉(zhuǎn)角導(dǎo)數(shù)相關(guān)法構(gòu)造出一個(gè)新的有效的厚薄板通用矩形單元，采用集中質(zhì)量矩陣法得到單元的質(zhì)量矩陣。通過(guò)編制MATLAB程序，采用子空間迭代法計(jì)算結(jié)構(gòu)在不同板厚時(shí)的自振頻率。將計(jì)算結(jié)果與理論解和大型通用有限元軟件ABAQUS解進(jìn)行比較，該單元具有良好的數(shù)值計(jì)算結(jié)果，能夠?qū)崿F(xiàn)動(dòng)力學(xué)范圍內(nèi)薄板、中厚板以及厚板完全通用，具有工程應(yīng)用價(jià)值。

關(guān)鍵詞：結(jié)構(gòu)力學(xué)；自振頻率；厚薄板通用單元；COONS曲面；廣義協(xié)調(diào)；MATLAB

對(duì)于結(jié)構(gòu)計(jì)算而言，靜力問(wèn)題分析雖是首要的，但是往往動(dòng)力作用引起的破壞才是致命的，這是引起結(jié)構(gòu)嚴(yán)重破壞的主因。例如：地震作用引起的建筑物的坍塌；風(fēng)振作用引起的高塔、橋梁的振動(dòng)破壞；碰撞或者是爆炸對(duì)結(jié)構(gòu)引起的沖擊破壞等等。因此，在工程結(jié)構(gòu)的研究、設(shè)計(jì)和安全性評(píng)估中，進(jìn)行結(jié)構(gòu)的振動(dòng)特性分析是相當(dāng)重要的，其中結(jié)構(gòu)振動(dòng)頻率是非常重要和基礎(chǔ)的動(dòng)力特性。當(dāng)前，板結(jié)構(gòu)是建筑結(jié)構(gòu)中最為常用也是最重要的結(jié)構(gòu)形式，因此對(duì)板結(jié)構(gòu)進(jìn)行自由振動(dòng)分析非常有必要。

由于板的特殊性，薄板與厚板的力學(xué)性能有所不同，對(duì)于薄板和厚板通常需分開(kāi)單獨(dú)計(jì)算，實(shí)際應(yīng)用中對(duì)于薄板與厚板的劃分無(wú)明顯界限。為了提高計(jì)算效率和精度，人們開(kāi)始探索厚薄板通用單元，而要實(shí)現(xiàn)厚薄板通用，則必須克服薄板情況下的剪切閉鎖問(wèn)題。此外，對(duì)于厚薄板殼的振動(dòng)問(wèn)題，采用常規(guī)的解析法求解是很困難的，而采用有限元法求解則能得到很好的效果[1-7]。有限元COONS曲面法作為有限單元法與計(jì)算幾何學(xué)相結(jié)合的產(chǎn)物，構(gòu)造單元簡(jiǎn)便高效，具有可以直接得到位移場(chǎng)函數(shù)顯式的優(yōu)點(diǎn)。

然而，到目前為止關(guān)于厚薄板通用的COONS曲面單元的振動(dòng)分析研究尚少。文獻(xiàn)[2—3]通過(guò)Timoshenko深梁理論推導(dǎo)剪切變形因子并從薄板理論出發(fā)構(gòu)造出廣義協(xié)調(diào)COONS曲面混合單元，能夠?qū)崿F(xiàn)動(dòng)力下的中厚薄板通用。由于基于薄板理論構(gòu)造，所以該單元缺乏轉(zhuǎn)角位移場(chǎng)和剪切應(yīng)變矩陣，對(duì)于厚跨比大的厚板(h/L>0.2時(shí))，計(jì)算結(jié)果并不理想。因此，本文從Mindlin-Reissner厚板理論出發(fā)并結(jié)合廣義協(xié)調(diào)思想，采用第2類(lèi)COONS曲面構(gòu)造法及撓度和轉(zhuǎn)角導(dǎo)數(shù)相關(guān)法構(gòu)造出一個(gè)新的厚薄板通用單元，通過(guò)編制MATLAB程序進(jìn)行振動(dòng)分析，為運(yùn)用COONS曲面法研究板殼振動(dòng)提供新的思路。

1單元構(gòu)造

彈性厚板理論(Mindlin-Reissner厚板理論)[8-9]假設(shè)原來(lái)垂直板中面的直線(xiàn)在變形后仍保持為直線(xiàn)，但由于考慮了剪切變形的影響，其不再垂直于變形后的板中面。因此，采用此理論的板單元撓度位移w和轉(zhuǎn)角位移θx，θy需各自獨(dú)立插值?；诤癜謇碚搧?lái)構(gòu)造厚薄板通用單元，關(guān)鍵在于如何在薄板情況下克服剪切閉鎖。為此，本文從單向厚板理論出發(fā)，推導(dǎo)出內(nèi)含考慮剪切變形影響的剪切變形因子λ的插值基函數(shù)。構(gòu)造單元時(shí)通過(guò)引入此插值基函數(shù)使板單元在薄板時(shí)的剪切剛度趨近于零，避免出現(xiàn)剪切閉鎖。

1.1剪切變形因子

如圖1所示，設(shè)單向厚板的撓度位移w為三次曲線(xiàn)，法線(xiàn)平均轉(zhuǎn)角θS為二次函數(shù)，則有式(1)：

(1)

受Timoshenko深梁理論推導(dǎo)的啟發(fā)，采用單向厚板理論推導(dǎo)得到[10-11]

(2)

圖1　單向厚板彎曲變形Fig.1　Bending deformation of one-way thick plate

將式(1)代入式(2)，積分得：

(3)

引入邊界條件w|S=0=wi,w|S=L=wj得：

C4=wi。

(4)

設(shè):

(5)

將式(4)和式(5)代入式(3)得單元撓度邊界

w=F0(t)wi+F1(t)wj+G0(t)θi+G1(t)θj。

(6)

將式(6)代入式(1)，得單元轉(zhuǎn)角邊界

θS=f0(t)wi+f1(t)wj+g0(t)θi+g1(t)θj。

(7)

式(6)、式(7)分別為單向厚板的撓度位移函數(shù)和轉(zhuǎn)角位移函數(shù)，均由物理方程推導(dǎo)得到。其中撓度位移插值基函數(shù)為

(8)

轉(zhuǎn)角位移插值基函數(shù)為

(9)

通過(guò)以上公式可以看出，當(dāng)板逐漸變薄時(shí)，即h/L→0時(shí)，λ→0，ρ→1時(shí)，

(10)

系數(shù)λ反映了剪切變形的影響，故稱(chēng)為剪切變形因子，它是薄板COONS曲面單元與厚板COONS曲面單元相互聯(lián)系的橋梁。

綜上所述，當(dāng)h/L→0時(shí)，

θS→w′。

(11)

1.2單元撓度位移場(chǎng)

進(jìn)行坐標(biāo)系轉(zhuǎn)換，無(wú)量綱自然坐標(biāo)系(ξ,η)與總體坐標(biāo)系(x,y)之間的關(guān)系為

ξ=x/a，η=y/b(0≤ξ,η≤1) ，

其中a和b為矩形單元的邊長(zhǎng)。

擬合單元4條邊界的撓度位移場(chǎng)，并保證其在撓度和切向?qū)?shù)方向上保持協(xié)調(diào)：

[w(ξ,0)w(ξ,1)]=

(12)

(13)

通過(guò)計(jì)算表明：以線(xiàn)性插值基函數(shù)為混合函數(shù)的ACM元的計(jì)算精度較高，但收斂速度較慢，而以第1類(lèi)Hermite插值基函數(shù)為混合函數(shù)的Ferguson曲面元收斂速度較快，但計(jì)算精度較差，據(jù)此本文對(duì)于單元邊界方向?qū)?shù)的擬合采用帶有剪切變形因子的線(xiàn)性組合函數(shù)為插值基函數(shù)進(jìn)行擬合，這樣可以保證該單元在收斂速度和計(jì)算精度上有所提高。擬合結(jié)果如下：

[wη(ξ,0)wη(ξ,1)]=

(14)

(15)

其中L0=1-t，L1=t，t=ξ,η。

再引入廣義協(xié)調(diào)條件，使COONS曲面片邊界處的法線(xiàn)導(dǎo)數(shù)滿(mǎn)足廣義協(xié)調(diào)，即令

(16)

[Nw]=

[F0F00.5b(L0+F0)G0-0.5aG0(L0+F0)

F1F00.5b(L1+F1)G0-0.5aG1(L0+F0)

F1F10.5b(L1+F1)G1-0.5aG1(L1+F1)

F0F10.5b(L0+F0)G1-0.5aG0(L1+F1)]，

{δ}=[w1θx1θy1w2θx2θy2

w3θx3θy3w4θx4θy4]T。

[w]=[Nw]{δ}。

(17)

[Nw]的各項(xiàng)乘積因式中，前者是ξ的函數(shù)，后者是η的函數(shù)，前后順序不可交換(下同)。

1.3單元轉(zhuǎn)角位移場(chǎng)

為避免剪切閉鎖，必須使厚板元的撓度位移場(chǎng)和轉(zhuǎn)角位移場(chǎng)在h/L→0時(shí)滿(mǎn)足導(dǎo)數(shù)關(guān)系，即滿(mǎn)足Kirchhoff薄板理論的直法線(xiàn)假設(shè)。式(10)和式(11)已指出，當(dāng)h/L→0時(shí)，撓度位移場(chǎng)w和轉(zhuǎn)角位移場(chǎng)θx，θy滿(mǎn)足導(dǎo)數(shù)關(guān)系。而h/L不趨于零時(shí)，撓度場(chǎng)和轉(zhuǎn)角場(chǎng)并不滿(mǎn)足導(dǎo)數(shù)關(guān)系，符合厚板理論的要求。因此，采用撓度位移場(chǎng)的“偏導(dǎo)數(shù)”作為轉(zhuǎn)角位移場(chǎng)，但求導(dǎo)時(shí)附加約束性求導(dǎo)條件[12]，

由此，就能得到厚薄板通用單元的轉(zhuǎn)角位移場(chǎng)函數(shù)：

具體表達(dá)式分別如式(18)和式(19)所示。

(18)

(19)

1.4單元?jiǎng)偠染仃?/p>

應(yīng)變有5個(gè)分量，即

{χ}=[χxχy2χxy…γxγy]T=

單元?jiǎng)偠染仃嚕?/p>

[k]e=[kb]e+[kS]e。

[kb]e=∫A[Bb]T[Db][Bb]dA為單元彎曲剛度矩陣，[kS]e=∫A[BS]T[DS][BS]dA為單元剪切剛度矩陣。

顯然，當(dāng)h/L→0時(shí)，[BS]→[0]，[kS]e→[0]，故不會(huì)出現(xiàn)剪切閉鎖。

1.5單元質(zhì)量矩陣

采用集中質(zhì)量矩陣：

其中ne為單元節(jié)點(diǎn)數(shù)，ρ為材料密度。

1.6編程求解

基于MATLAB[13]強(qiáng)大的矩陣運(yùn)算能力，本文采用MATLAB編制了結(jié)構(gòu)計(jì)算程序，并采用子空間迭代法求解結(jié)構(gòu)的自振頻率。計(jì)算中只需修改板厚，即可實(shí)現(xiàn)厚薄板的自振分析，無(wú)須像ABAQUS分析時(shí)需根據(jù)厚跨比選用不同的單元進(jìn)行計(jì)算。

2數(shù)值算例

數(shù)值算例描述：四邊簡(jiǎn)支各向同性矩形彈性方板，邊長(zhǎng)為10 m，材料密度ρ=2 400 kg/m3，彈性模量E=3.0×1010N/m2，泊松比v=1/6。采用MATLAB編制的程序，按12×12劃分單元，求解不同厚跨比(h/L)下前5階頻率值，并與理論解和ABAQUS解[14]進(jìn)行對(duì)照，結(jié)果如表1和圖2—圖6所示。

圖2　不同厚跨比一階頻率比較Fig.2　Comparison of the first order frequency 　　　by different thickness span ratio

圖3　不同厚跨比二階頻率比較Fig.3　Comparison of the second order frequency　　　by different thickness span ratio

階次頻率/Hz薄板h/L=0.001薄板h/L=0.01中厚板h/L=0.1中厚板h/L=0.15中厚板h/L=0.2中厚板h/L=0.25厚板h/L=0.3厚板h/L=0.35本文解0.32493.247430.70143.75655.35065.61074.66182.626(1,1)理論解0.32523.252031.44845.56557.96368.64077.59785.192ABAQUS解0.32713.268330.03742.71053.76763.33371.60278.763本文解0.81158.110075.116104.530128.290147.120161.930171.560(1,2)理論解0.81308.130175.451104.520126.860143.750156.180165.440ABAQUS解0.83348.326173.473101.450123.620141.030154.720165.540本文解1.295012.9370114.83154.520184.160206.180222.620234.980(2,2)理論解1.300813.0080116.04155.350182.970201.740214.800223.770ABAQUS解1.330912.2890111.590149.320177.230198.060204.070207.700本文解1.622216.2020143.500191.480225.190248.610265.040273.790(1,3)理論解1.626116.260141.320185.910215.390234.720247.550256.970ABAQUS解1.738517.353140.200186.430207.980210.850229.580231.740本文解2.100620.968177.750230.750266.200290.070306.490318.080(2,3)理論解2.114821.147177.620227.780259.350278.790291.450300.110ABAQUS解2.220122.150173.260208.730219.160245.850265.140278.400

注：表1中薄板的理論解參考文獻(xiàn)[15]，中厚板及厚板的理論解參考文獻(xiàn)[10]；ABAQUS模擬計(jì)算中，薄板采用S4R5單元，中厚板和厚板采用S8R單元。

圖4　不同厚跨比三階頻率比較Fig.4　Comparison of the three order frequency 　　　by different thickness span ratio

圖5　不同厚跨比四階頻率比較Fig.5　Comparison of the four order frequency　　　by different thickness span ratio

圖6　不同厚跨比五階頻率比較Fig.6　Comparison of the five order frequency　　　by different thickness span ratio

由以上表1和圖2—圖6可以看出：本單元在不同厚跨比下的前3階數(shù)值解基本與理論解吻合。隨著階數(shù)和板厚的增加，本文解的精度仍?xún)?yōu)于ABAQUS解，能夠滿(mǎn)足工程計(jì)算需要。本單元能夠在動(dòng)力學(xué)范圍內(nèi)實(shí)現(xiàn)薄板-中厚板-厚板的完全通用，能夠得到收斂速度快、計(jì)算精度較高的數(shù)值結(jié)果。

隨板厚和階數(shù)增加而出現(xiàn)的誤差與單元推導(dǎo)及程序編制時(shí)未考慮轉(zhuǎn)動(dòng)慣量和擠壓變形的影響有關(guān)系。但相比于剪切變形的影響，以上的影響很小。為了提高精度，以后可在這方面進(jìn)行相關(guān)研究。

3結(jié)論

基于Mindlin-Reissner厚板理論并結(jié)合廣義協(xié)調(diào)思想，采用第2類(lèi)COONS曲面構(gòu)造法及撓度和轉(zhuǎn)角導(dǎo)數(shù)相關(guān)法構(gòu)造出一個(gè)新的厚薄板通用單元，通過(guò)自編MATLAB程序進(jìn)行振動(dòng)分析，數(shù)值算例表明該單元具有良好的收斂性和計(jì)算精度，能夠?qū)崿F(xiàn)動(dòng)力學(xué)范圍內(nèi)薄板、中厚板及厚板的完全通用，具有工程實(shí)用性。

1)選用單向厚板基本方程推導(dǎo)出剪切變形因子，基于厚板理論構(gòu)造出的COONS曲面單元能夠?qū)崿F(xiàn)動(dòng)力下的完全厚薄板通用，可作進(jìn)一步的研究和推廣，以豐富大型有限元分析軟件的單元庫(kù)。

2)相比于第1類(lèi)COONS曲面，第2類(lèi)COONS曲面較復(fù)雜，自由度數(shù)也會(huì)增加，但收斂速度更快，計(jì)算精度更高。本文采用第2類(lèi)COONS曲面構(gòu)造法，通過(guò)引入廣義協(xié)調(diào)條件消去自由度，從而得到消去節(jié)點(diǎn)扭矢的12參數(shù)矩形元，同時(shí)在單元邊界法向?qū)?shù)擬合中使用內(nèi)含剪切變形因子的線(xiàn)性組合函數(shù)作為插值基函數(shù)。數(shù)值算例結(jié)果表明，本單元具有良好的數(shù)值計(jì)算結(jié)果。

3)COONS曲面單元具有直接給出位移場(chǎng)顯式的優(yōu)點(diǎn)，已在靜力問(wèn)題分析中取得了較好的收斂性和精度，然而在動(dòng)力問(wèn)題分析中缺少研究。本文通過(guò)編制的MATLAB程序，對(duì)推導(dǎo)出的基于厚板理論的厚薄板通用的COONS曲面單元進(jìn)行振動(dòng)特性分析。計(jì)算表明，COONS曲面單元在動(dòng)力問(wèn)題分析中同樣能夠?qū)崿F(xiàn)厚薄板通用，且具有良好的數(shù)值計(jì)算結(jié)果，為運(yùn)用有限元COONS曲面法進(jìn)行板殼的振動(dòng)分析提供了借鑒。

參考文獻(xiàn)/References：

[1]張延慶.板殼有限元計(jì)算幾何理論[D]. 徐州:中國(guó)礦業(yè)大學(xué), 1990.

ZHANG Yanqing. Computational Geometric Theory of Finite Element of Plate and Shell[D]. Xuzhou: China University of Mining and Technology，1990.

[2]張延慶，龍馭球. 厚薄板通用的廣義協(xié)調(diào)單元COONS曲面構(gòu)造法[J]. 工程力學(xué), 1994, 11(2)：59-64.

ZHANG Yanqing, LONG Yuqiu. Formulation of a generalized conforming element for moderate thick plate based on COONS’ surface[J]. Engineering Mechanics，1994，11(2)：59-64.

[3]陳葉妮，張延慶. 用廣義協(xié)調(diào)COONS混合曲面元分析中厚板的振動(dòng)[J]. 工業(yè)建筑, 2006(sup1)：612-614.

CHEN Yeni，ZHANG Yanqing. Vibration analysis of moderate thick plate by generalized conforming element and COONS surface approach[J]. Industrial Construction，2006(sup1)：612-614.

[4]曹金鳳，姜耀東，趙國(guó)景，等. 廣義協(xié)調(diào)厚薄板通用單元TMT的固有振動(dòng)分析[J]. 中國(guó)礦業(yè)大學(xué)學(xué)報(bào)，2007, 36(1)：91-96.

CAO Jinfeng，JIANG Yaodong，ZHAO Guojing，et al. Free vibration analysis of plate using a generalized conforming thick/thin element TMT[J]. Journal of China University of Mining and Technology，2007, 36(1)：91-96.

[5]曹金鳳，姜耀東，韓志茹，等. 廣義協(xié)調(diào)厚薄板通用單元的固有振動(dòng)分析[J]. 力學(xué)與實(shí)踐，2007，29(3)：30-35.

CAO Jinfeng，JIANG Yaodong，HAN Zhiru，et al. Free vibration analysis for plates using a generalized conforming thick/thin triangular element[J]. Journal of Mechanics in Engineering，2007, 29(3)：30-35.

[6]韓志茹，姜耀東，王麗，等. 廣義協(xié)調(diào)厚薄通用板單元的動(dòng)力特性分析[J]. 力學(xué)與實(shí)踐，2009, 31(3)：69-72.

HAN Zhiru，JIANG Yaodong，WANG Li，et al. The dynamic property of generalized conforming thick/thin plate elements[J]. Journal of Mechanics in Engineering，2009, 31(3)：69-72.

[7]黃修長(zhǎng)，錢(qián)振華，李俊，等. 應(yīng)用基于解析式函數(shù)的廣義協(xié)調(diào)四邊形厚板元分析中厚板的自由振動(dòng)[J]. 工程力學(xué)，2011, 28(9)：39-43.

HUANG Xiuchang，QIAN Zhenhua，LI Jun，et al. Vibration analysis of moderate thick plate by generalized conforming quadrilateral thick plate element based on analytical trial functions[J]. Engineering Mechanics，2011, 28(9)：39-43.

[8]MINDLIN R D. Influence of rotatory inertia and shear on flexural motions of isotropic elastic plates[J]. ASME Journal of Applied Mechanics，1951, 18(1)：31-38.

[9]REISSNER E. On bending of elastic plates[J]. Quarterly of Applied Mathematics，1947, 5(1)：55-68.

[10]曹志遠(yuǎn)，楊昇田. 厚板動(dòng)力學(xué)理論及其應(yīng)用[M]. 北京:科學(xué)出版社，1983：57-124.

[11]鄭照北，王宏凱. 板殼有限元COONS曲面法[M]. 徐州:中國(guó)礦業(yè)大學(xué)出版社，1993：52-73.

[12]王宏凱，鄭照北. 一個(gè)有效的厚薄板通用COONS曲面矩形元[J]. 工程力學(xué)，1992, 9(1)：134-140.

WANG Hongkai，ZHENG Zhaobei. An efficacious COONS’ surface method for formulation of current finite element of thick and thin plates[J]. Engineering Mechanics，1992, 9(1)：134-140.

[13]劉衛(wèi)國(guó). MATLAB程序設(shè)計(jì)教程[M].第2版. 北京：中國(guó)水利水電出版社，2010.

[14]莊茁，曲小川，廖劍暉, 等. 基于ABAQUS的有限元分析和應(yīng)用[M]. 北京：清華大學(xué)出版社，2009.

[15]曹志遠(yuǎn). 板殼振動(dòng)理論[M]. 北京: 中國(guó)鐵道出版社, 1989：32-54.

Vibration analysis of generalized conforming and COONS’surface thick/thin plate element

CHEN Bin, ZHANG Yanqing

(College of Architecture and Civil Engineering, Beijing University of Technology, Beijing 100124, China)

Abstract:Thick/thin plate element has been the focus and difficulty in the study of shell finite element. According to the basic equation of one-way thick plate, the base function of interpolation containing the shear deformation factor is derived. Based on the Mindlin-Reissner thick plate theory and the idea of generalized coordination, the second kind of finite element COONS’ surface construction method and the derivative correlation method of deflection displacement and angular rotation displacement are used to derive a new effective thick/thin general rectangular plate element, then the lumped mass matrix method is used to get the mass matrix of the element. By using MATLAB programming, the natural frequencies of the plates with different thickness are acquired by a sub-space iteration method. Compared with the theoretical solution and the numerical value of ABAQUS, it is concluded that this element has good numerical calculation results, which can realize entirely common use for thin plate, cut deal and thick plate in dynamic range, and it has engineering application value.

Keywords:structural mechanics; natural frequency; thick/thin plate element; COONS’ surface; generalized conforming; MATLAB

文章編號(hào)：1008-1534(2016)03-0224-06

收稿日期：2016-04-04;修回日期：2016-04-24；責(zé)任編輯：馮民

作者簡(jiǎn)介：陳斌(1990—)，男，浙江金華人，碩士研究生，主要從事工程結(jié)構(gòu)數(shù)值分析與優(yōu)化方面的研究。通訊作者：張延慶教授。E-mail:zhyq@bjut.edu.cn

中圖分類(lèi)號(hào)：TU311.4

文獻(xiàn)標(biāo)志碼：A

doi:10.7535/hbgykj.2016yx03008

陳斌，張延慶.廣義協(xié)調(diào)COONS曲面厚薄板通用單元振動(dòng)分析[J].河北工業(yè)科技,2016,33(3):224-229.

CHEN Bin, ZHANG Yanqing.Vibration analysis of generalized conforming and COONS’ surface thick/thin plate element[J].Hebei Journal of Industrial Science and Technology,2016,33(3):224-229.