混合能量原理求解矩形板在靜水壓力作用下的彎曲

陳英杰 闕春發(fā) 郭敦

(燕山大學(xué)建筑工程與力學(xué)學(xué)院,河北秦皇島066004)

矩形板作為結(jié)構(gòu)基本構(gòu)件之一，我們可以把矩形板理論分析方法歸納分為精確解法和近似解法，且這兩類方法均以能量原理為基礎(chǔ)。能量法又稱為變分法，可以追溯到17世紀(jì)末，有著悠久的歷史。變分法可以用于處理函數(shù)變量的問題。截止到目前，國內(nèi)外學(xué)者普遍非常重視變分原理[1-2]的研究，而且已經(jīng)有了豐碩的成果。例如：王根會(huì)等[3]基于能量變分法推導(dǎo)了新型組合箱梁的總勢能式與微分控制方程，為新型組合箱梁的應(yīng)用推廣提供了理論和技術(shù)支持；曾祥勇等[4]以能量變分法的最小勢能原理為基礎(chǔ)分析了矩形筏板的受彎問題?；旌献兞康臉O值變分原理[5-6]包含混合變量的最小勢能原理[7-9]與混合變量的最小余能原理[10]。混合變量的極值變分原理，其容許位移與容許內(nèi)力分別被弱容許位移與弱容許內(nèi)力所代替，與傳統(tǒng)理論[11-12]相比增強(qiáng)了等價(jià)方程，而且需要預(yù)先滿足的條件減弱了。

本文將應(yīng)用變分原理中的混合變量的最小勢能原理討論矩形板在靜水壓力作用下的彎曲問題，并最終得到矩形板受靜水壓力的數(shù)值精確解。應(yīng)用Matlab軟件計(jì)算出矩形板撓度的精確值，應(yīng)用ANSYS模擬軟件得到矩形板撓度模擬值，并對(duì)兩項(xiàng)數(shù)值進(jìn)行對(duì)比分析，說明本文方法具有更好的準(zhǔn)確性。

1 三邊簡支一邊固定矩形板

1.1 撓曲線方程

現(xiàn)考慮受靜水壓力的矩形板，邊界條件如圖1(a)所示。圖中q0為三角形分布的均布載荷，b為矩形板沿y軸方向的長度，a為矩形板沿x軸方向的長度，w為撓度。以分布彎矩My0代替固定邊的彎曲約束，創(chuàng)建如圖1(b)所示的變形后的等效圖。

圖1 受靜水壓力的三邊簡支一邊固定的矩形板

分布載荷

假設(shè)固定端彎矩

式中Cm為待定系數(shù)，αm=mπ/a，同時(shí)假設(shè)板的弱容許撓度為

式中，βn=nπ/b。

圖1(b)對(duì)應(yīng)的矩形板混合變量總勢能

式中D=Eh3/[12(1?ν2)]，為板的抗彎剛度，E為彈性模量，ν為泊松比。

將式(1)和式(2)代入式(4)經(jīng)過積分運(yùn)算并且對(duì)Amn取變分極值得

將式(5)代入式(3)得

在推導(dǎo)本邊界條件撓曲線方程的過程中，式(6)是用混合變量最小勢能原理導(dǎo)出的以正弦雙重三角級(jí)數(shù)表示的弱容許撓度。而對(duì)于非齊次撓度和彎矩邊界條件，正弦雙重三角級(jí)數(shù)將會(huì)在邊界上出現(xiàn)第二類間斷點(diǎn)。為避免出現(xiàn)這種情況，并且能夠加快級(jí)數(shù)收斂速度，需要將式(6)轉(zhuǎn)換成為在邊界上連續(xù)可微的撓度，而該撓度即為本次邊界條件下的撓度表達(dá)式。應(yīng)用文獻(xiàn)[5]中附錄式(A92)和式(A47)分別對(duì)式(6)各項(xiàng)進(jìn)行轉(zhuǎn)化，則得

1.2 邊界條件

固定端轉(zhuǎn)角為0，故式(7)應(yīng)符合式(8)的邊界條件，其他邊界條件已自動(dòng)滿足

將方程式(7)代入式(8)中，經(jīng)過運(yùn)算整理可得執(zhí)行方程

1.3 數(shù)值計(jì)算

取矩形板的各項(xiàng)參數(shù)：a=b=1 m，泊松比ν=0.3,彈性模量E=2.0×1011Pa，板的厚度h=0.01 m，靜水壓力q0=1 MPa。其中含有未知項(xiàng)Cm。其次，應(yīng)用Matlab軟件編程計(jì)算，可求得具體撓度值。對(duì)計(jì)算程序循環(huán)50次，便可保障計(jì)算結(jié)果的收斂性。

同時(shí)，應(yīng)用模擬軟件對(duì)本次邊界條件下的矩形板進(jìn)行建模分析，作為參考解與本文解進(jìn)行對(duì)比。此種邊界條件選用模型shell63單元，該單元具有彎曲能力和薄膜效應(yīng)，并且忽略剪切變形，可以很好地應(yīng)用于矩形板變形過程中的數(shù)值模擬。彎曲矩形板在x/a=0.1，x/a=0.3，x/a=0.5和x/a=0.8與y/b=0.1，y/b=0.3，y/b=0.5和y/b=0.8線上的模擬數(shù)值與本文計(jì)算數(shù)值兩項(xiàng)數(shù)據(jù)繪制成圖和表，如圖2與表1所示。

圖2 x/a=0.1,x/a=0.3,x/a=0.5和x/a=0.8處撓度分布曲線圖

1.4 結(jié)果分析

圖2為矩形板沿x軸不同位置處的本文計(jì)算值和模擬值的撓度對(duì)比曲線圖，可以直觀看到,由于結(jié)構(gòu)受三角形載荷作用，矩形板在靜水壓力作用下的撓度并非對(duì)稱分布變化，且符合實(shí)際受力變化規(guī)律，計(jì)算結(jié)果更為貼近實(shí)際，這也間接說明混合變量的最小勢能原理解決矩形板受靜水壓力問題的適用性。

表1給出了矩形板在y/b=0.1，y/b=0.3，y/b=0.5和y/b=0.8位置處沿x軸方向變化的撓度計(jì)算值和模擬參考值。通過對(duì)兩項(xiàng)數(shù)值的對(duì)比分析，可以得到兩者撓度的最大相對(duì)差值分別為：4.0%，3.8%，3.7%，本文研究方法計(jì)算的撓度值略大于有限元模擬的參考解，表明本方法可以更好地保證結(jié)構(gòu)安全，反映結(jié)構(gòu)在載荷作用下的變形規(guī)律，對(duì)矩形板彎曲解擁有更好的有效性和適用性，可以更精準(zhǔn)地計(jì)算矩形板的撓度問題。

表1 y/b=0.1，y/b=0.3，y/b=0.5和y/b=0.8處撓度沿x分布值(10?7 m)

2 兩鄰邊固定兩鄰邊簡支矩形板

2.1 撓曲線方程

現(xiàn)考慮受靜水壓力的矩形板，邊界條件如圖3(a)所示。以分布彎矩Mx0，My0代替固定邊的彎曲約束，創(chuàng)建如圖3(b)所示的變形后的等效圖。

圖3 受靜水壓力的兩鄰邊固定兩鄰邊簡支矩形板

分布載荷

假設(shè)固定端彎矩

式中，An為待定系數(shù)。同時(shí)假設(shè)板的弱容許撓度為

可得圖3(b)對(duì)應(yīng)的矩形板混合變量總勢能式

將式(10)～式(12)代入式(14)經(jīng)過積分運(yùn)算并且對(duì)Amn取變分極值得

將式(15)代入式(13)得

同理，需要將式(16)轉(zhuǎn)換成雙曲函數(shù)和三角級(jí)數(shù)混合表示的基本解。應(yīng)用文獻(xiàn)[5]中附錄式(A47)和式(A92)分別對(duì)式(16)各項(xiàng)進(jìn)行轉(zhuǎn)化，得

2.2 邊界條件

固定端應(yīng)滿足轉(zhuǎn)角為零，故式(17)應(yīng)符合式(18)和式(19)的邊界條件，其他邊界條件已自動(dòng)滿足

將方程式(17)代入式(18)和式(19)中，運(yùn)算整理得執(zhí)行方程

2.3 數(shù)值計(jì)算

通過之前的公式推導(dǎo)，已經(jīng)得到本邊界條件的執(zhí)行方程。且矩形板的計(jì)算參數(shù)與1.3相同。本節(jié)中所需要求解的未知數(shù)是An和Cm，對(duì)計(jì)算程序循環(huán)50次，運(yùn)行結(jié)果不再變化，得到本邊界條件下的撓度值。應(yīng)用模擬軟件對(duì)本次邊界條件下的矩形板進(jìn)行建模分析，單元類型采用Shell63，將彎曲矩形板在x/a=0.1，x/a=0.3，x/a=0.5和x/a=0.8與y/b=0.1，y/b=0.3，y/b=0.5和y/b=0.8線上的模擬數(shù)值與本文計(jì)算數(shù)值分別繪制成圖和表，如圖4與表2所示。

圖4 x/a=0.1，x/a=0.3，x/a=0.5和x/a=0.8處撓度分布曲線圖

2.4 結(jié)果分析

圖4給出了矩形板沿x軸不同位置處的本文計(jì)算值和模擬值的撓度對(duì)比曲線圖，可以直觀看出，本文研究方法計(jì)算的撓度值略大于有限元模擬的參考解，表明本方法可以更好地保證結(jié)構(gòu)安全，反映結(jié)構(gòu)在載荷作用下的變形規(guī)律，對(duì)矩形板彎曲解具有更好的有效性和適用性，可以更精準(zhǔn)地計(jì)算矩形板的撓度問題。由于受到載荷分布位置與邊界條件不同的影響，矩形板的撓度曲線均沿坐標(biāo)軸先增大后減小，且在x/a=0.6的位置達(dá)到了最大值。表明計(jì)算出的矩形板在靜水壓力作用下的撓度變化是符合受力變化規(guī)律的。

表2給出了矩形板在y/b=0.1，y/b=0.3，y/b=0.5和y/b=0.8位置處沿x軸方向變化的撓度公式計(jì)算值和模擬值。通過對(duì)兩項(xiàng)數(shù)值的對(duì)比分析，可以得到兩者撓度的最大相對(duì)差值分別為：3.5%、3.4%、3.0%，本文研究方法計(jì)算的撓度值略大于有限元模擬的參考解，表明本方法可以更好地保證結(jié)構(gòu)安全，反映結(jié)構(gòu)在載荷作用下的變形規(guī)律，使矩形板撓度問題的計(jì)算更為精確。

表2 矩形板在y/b=0.1，y/b=0.3，y/b=0.5和y/b=0.8處撓度沿x分布值(10?7 m)

3 三邊固定一邊簡支矩形板

3.1 撓曲線方程

現(xiàn)考慮受靜水壓力的矩形板，邊界條件如圖5(a)所示。以分布彎矩Mx0，My0，Mxa代替固定邊的彎曲約束，創(chuàng)建如圖5(b)所示的變形后的等效圖。

分布載荷

假設(shè)固定端彎矩

假設(shè)板的弱容許撓度為

可得圖5(b)對(duì)應(yīng)的矩形板混合變量總勢能式

圖5 受靜水壓力的三邊固定一邊簡支彎曲矩形板

將式(22)～式(24)代入式(26)經(jīng)過積分運(yùn)算并且對(duì)Amn取變分極值得

將式(27)代入式(25)得

同理，需要將式(28)轉(zhuǎn)換成雙曲函數(shù)和三角級(jí)數(shù)混合表示的基本解。應(yīng)用文獻(xiàn)[5]中附錄式(A47)和式(A92)分別對(duì)式(28)各項(xiàng)進(jìn)行轉(zhuǎn)化，可得

3.2 邊界條件

固定端應(yīng)滿足轉(zhuǎn)角為零，故式(29)應(yīng)符合式(30)和式(31)的邊界條件，其他邊界條件已自動(dòng)滿足

把方程式(29)代入式(30)和式(31)中，經(jīng)過運(yùn)算整理得執(zhí)行方程。

3.3 數(shù)值計(jì)算

通過之前的公式推導(dǎo)，已經(jīng)得到本邊界條件的執(zhí)行方程。且矩形板的計(jì)算參數(shù)與1.3相同。本節(jié)中所需要求解的未知數(shù)是An和Cm，同理為了保障計(jì)算數(shù)值的收斂性，對(duì)計(jì)算程序循環(huán)50次，得到本邊界條件下的撓度值，并通過模擬軟件模擬分析，將模擬結(jié)果作為參考值。將彎曲矩形板在x/a=0.1，x/a=0.3，x/a=0.5和x/a=0.8與y/b=0.1，y/b=0.3，y/b=0.5和y/b=0.8線上的模擬數(shù)值與本文計(jì)算數(shù)值兩項(xiàng)數(shù)據(jù)分別繪制成圖和表，如圖6與表3所示。

3.4 結(jié)果分析

圖6給出了矩形板沿y軸不同位置處的本文計(jì)算值和模擬值的撓度對(duì)比曲線圖，可以直觀看出，由于受到載荷分布位置與邊界條件不同的影響，矩形板的撓度曲線都沿坐標(biāo)軸先增大后減小，且在x/a=0.6的位置達(dá)到了最大值。表明計(jì)算出的矩形板在靜水壓力作用下的撓度變化符合受力變化規(guī)律。

圖6 x/a=0.1，x/a=0.3，x/a=0.5和x/a=0.8處撓度分布曲線圖

表3給出了矩形板在y/b=0.1，y/b=0.3，y/b=0.5和y/b=0.8位置處沿x軸方向變化的撓度公式計(jì)算值和模擬值。通過對(duì)兩項(xiàng)數(shù)值的對(duì)比分析，可以得到兩者撓度的最大相對(duì)差值分別為：4.2%、3.9%、3.2%，本文研究方法計(jì)算的撓度值略大于有限元模擬的參考解，表明本方法可以更好地保證結(jié)構(gòu)安全，并且對(duì)求解矩形板彎曲解擁有更好的有效性和適用性。

表3 y/b=0.1，y/b=0.3，y/b=0.5和y/b=0.8處撓度沿x分布值(10?7 m)

4 結(jié)論

本文依據(jù)混合變量的最小勢能原理推導(dǎo)出了三邊簡支一邊固定、兩鄰邊固定兩鄰邊簡支、三邊固定一邊簡支三種不同邊界條件的彎曲矩形板在靜水壓力作用下的邊界應(yīng)力函數(shù)的表達(dá)式和封閉解析解，數(shù)值計(jì)算軟件求得的數(shù)值解與模擬軟件的模擬值進(jìn)行歸納分析，表明本文方法的正確性，說明了本文研究方法可以使矩形板的彎曲問題的求解簡單化和精確化。