求解加筋板的位移與應(yīng)力的一種解析計算方法

王小明

(中國艦船研究設(shè)計中心 武漢 430064)

0 引 言

加筋板是船舶、鋼結(jié)構(gòu)橋梁和海洋平臺的主要構(gòu)件.在橋梁領(lǐng)域，由于其縱向和橫向加筋導(dǎo)致這兩個方向的等效剛度不相等，通常將縱橫加筋板稱為正交異性板.對于這類正交異性板的計算理論有三個體系：①主梁體系，由面板、縱筋和橫筋組成主向梁，抵抗整體變形；②梁板體系，由面板、縱筋和橫筋組成一個整體，面板作為縱筋和橫筋的共同翼板；③面板體系，面板視為支撐在縱筋和橫筋上的各向同性連續(xù)板[1-3].對于第一體系，確定了翼板寬度，其力學(xué)分析問題可以由初等梁理論解決；對于第三體系，可以采用薄板理論分析；對于第二體系，其存在梁板耦合的系統(tǒng)，計算比較復(fù)雜，正是本文研究的內(nèi)容.

關(guān)于第二體系變形和應(yīng)力的計算方法，Younseok等[4]提出了將加筋板等效成正交異性板的詳細(xì)計算方法.這樣的算法存在一定的局限性，要求加強(qiáng)筋比較小，且密集布置，載荷均布或者加載在加強(qiáng)筋位置上，不滿足這兩點(diǎn)，計算精度就比較差，并且無法計算加強(qiáng)筋上的應(yīng)力和位移.Jeom等[5]提出了基于Huber理論的P-E法，將縱橫加筋板的彎曲計算分為兩個步驟：①假設(shè)橫向加強(qiáng)筋剛度無限大，按剛性支持的正交異性板計算；②計算橫向加強(qiáng)筋的影響，然后修正第一步的計算結(jié)果.Bhaskar等[6]采用經(jīng)典薄板理論和歐拉-伯努利梁理論研究了各向同性加筋板和正交異性加筋板的彎曲和自由振動，該方法的計算結(jié)果與三維有限元結(jié)果吻合良好.采用一階剪切變形理論和能量變分法推導(dǎo)加筋板的自由振動方程是常用的方法，關(guān)于振動微分方程的解法，Naveed等[7]并沒有完全采用廣為熟悉的雙三角級數(shù)的方法，而是引入了鐵木辛科(Timoshenko)梁函數(shù)作為形函數(shù)；Hamedani等[8]則是把微分方程轉(zhuǎn)化為有限元格式求解.Shanmugam等[9]采用有限元計算和試驗驗證的方式研究了加筋板在面內(nèi)載荷和橫向載荷作用下的變形和破壞，研究結(jié)果認(rèn)為面內(nèi)載荷增加，橫向承載力將下降.彭林欣等[10-11]基于移動最小二乘近似，給出了加肋板線性彎曲的無網(wǎng)格法.運(yùn)用無網(wǎng)格解得出的結(jié)果與加肋板三點(diǎn)彎曲實驗的位移結(jié)果吻合良好，與ANSYS有限元求解結(jié)果誤差不超過5%.Sapountzakis等[12-14]提出了一種板筋隔離法計算加筋板彎曲問題，將板筋拆解后，獨(dú)立研究，認(rèn)為板筋系統(tǒng)的內(nèi)力連同真實外力看作板的外力，板筋內(nèi)力看作筋的外力，應(yīng)用板筋結(jié)合面上的位移連續(xù)條件，求解橫向位移.研究中解微分方程的方法轉(zhuǎn)化為有限元格式計算，算例中重點(diǎn)關(guān)注的是橫向位移和板筋內(nèi)力分布，對應(yīng)力計算結(jié)果沒有研究.Lin等[15-17]先后用經(jīng)典薄板理論、歐拉-伯努利梁理論、一階剪切變形理論和鐵木辛科梁理論推導(dǎo)加筋板的振動微分方程，以梁板間的耦合力和力矩為中間未知量，采用位移協(xié)調(diào)條件求解了自由邊界條件下、固支、簡支邊界條件下板的速度導(dǎo)納.在求解中，沒有考慮到筋板之間的水平剪切力，認(rèn)為板筋的彎曲中面都是自身的幾何中面，在加強(qiáng)筋尺寸比較小，分布稀疏時計算結(jié)果差別不大，若加強(qiáng)筋尺寸大，且分布密集，則誤差較大.其重點(diǎn)討論的是一根筋的情況，實際工程中，往往是多根筋的，因此板筋之間的耦合力和力矩在不同加強(qiáng)筋之間互相影響；研究中關(guān)于加筋板的位移應(yīng)力并沒有討論.雖然有限元計算已經(jīng)成為當(dāng)前的主流算法，但是采用通用有限元軟件計算，計算結(jié)果與網(wǎng)格劃分、計算者的水平和經(jīng)驗密切相關(guān).工程中的設(shè)計過程，經(jīng)常是方案不斷變更，方案變更導(dǎo)致有限元模型修改工作量巨大.所以，在船舶設(shè)計領(lǐng)域，大部分的強(qiáng)度和剛度計算依然是以經(jīng)典理論公式計算為主.在這樣的背景下，研究加筋板的解析計算方法依然非常必要.

為了求解加筋板的彎曲問題，文中提出了一種針對多根加強(qiáng)筋構(gòu)成的加筋板的“板筋分離法”，與文獻(xiàn)[12-14]不同的是所有求解過程均采用半解析算法.這種方法以板筋之間的作用力和力矩為中間未知變量，將板筋之間的力、力矩和真實外力看作加筋板中板的外力，將板筋之間的力和力矩看作筋的外力，分別獨(dú)立求解二者的彎曲位移，通過位移協(xié)調(diào)條件確定板筋之間的耦合力和力矩.最終確定加筋板的位移和應(yīng)力.這種方法實質(zhì)上考慮了板筋之間的剪切力和法向力，比采用正交異性板法精度高；避免了能量法中需要針對不同截面分別計算變形能的麻煩，因此比能量法具有更廣的應(yīng)用范圍，可以應(yīng)用于工字型鋼加筋板，倒T型鋼加筋板等.

1 基本微分方程

建立坐標(biāo)系，見圖1.加筋板沿x軸方向長度為a，沿y軸方向?qū)挾葹閎，板厚為t.加強(qiáng)筋為角鋼，厚度為ts，高度為h.加強(qiáng)筋平行于x軸方向布置.外力大小為F0，作用點(diǎn)的坐標(biāo)為(x0，y0).加筋板橫向彎曲時，自由度可以用φx，φy和w表示.φx為xz平面內(nèi)板的轉(zhuǎn)角，φy為yz平面內(nèi)板的轉(zhuǎn)角，都以板中面轉(zhuǎn)向z軸的方向為正.w為板的垂向(z向)位移，u0p，v0p分別為板中面的水平(x向和y向)位移.板和筋是焊接在一起的，因此分界面上具有相同的變形.

圖1 加筋板與坐標(biāo)系統(tǒng)

受文獻(xiàn)[12-13]的啟發(fā)，并考慮板筋之間剪切力引起的剛體位移，對于多根加強(qiáng)筋的情形，板的彎曲微分方程可表示為

(1)

(2)

(3)

(4)

(5)

式中：N為加強(qiáng)筋的總數(shù)；D=Et3/12(1-μ2)為板的抗彎剛度；C=κEt/2(1+μ)為抗剪剛度；κ=5/6為剪切修正系數(shù)；E，μ分別為板的彈性模量和泊松比；δ為狄拉克函數(shù).

四邊簡支的加筋板，其邊界條件可以表示為

x=0或a，則

(6)

y=0或b，則

(7)

對于加強(qiáng)筋而言，將其看作鐵木辛科梁，其自由度為xz平面的彎曲轉(zhuǎn)角θ(x)，yz平面的扭轉(zhuǎn)轉(zhuǎn)角φs(x)和垂向(z向)位移u(x).u0s為梁中面的水平(x向)位移.

以加強(qiáng)筋形心點(diǎn)為原點(diǎn)建立坐標(biāo)系，即C為角鋼橫截面的形心位置，見圖2.加筋板在橫向彎曲時，第i根加強(qiáng)筋受到板對它的力qi、qxi和轉(zhuǎn)矩Mi的作用.

圖2 加筋板中板梁的耦合力與力矩

圖2中的Mi，僅示意其大小分布，方向由板(筋)上的雙箭頭表示.所以第i根加強(qiáng)筋的控制微分方程為

(8)

(9)

(10)

(11)

簡支梁的邊界條件為

x=0或a，則

(12)

2 基本微分方程的解

由式(1)～(5)組成的方程組與邊界條件式(6)～(7)構(gòu)成定解問題.可以采用雙傅里葉級數(shù)來求解這一定解問題.

(16)

(17)

式中：km=mπ/a,kn=nπ/b.式(13)～(17)已經(jīng)滿足了邊界條件式(6)～(7)，則僅需要通過方程組確定傅里葉級數(shù)的系數(shù).方程(1)～(3)與方程(4)～(5)是互相獨(dú)立的，可以分別求解.將式(13)～(15)代入方程式(1)～(3).

(18)

(19)

由式(18)可以分別求出傅里葉系數(shù).

(20)

(21)

(22)

方程(4)和(5)也可以采用相同的方法計算.

(23)

(24)

(25)

(26)

對加強(qiáng)筋所在的梁也可以作類似的求解與推導(dǎo).以第i根加強(qiáng)筋為例，方程式(8)～(11)與邊界條件式(12)組成的定解問題，解的形式可以寫成

(27)

(28)

(29)

(30)

式(27)～(30)代入式(8)～(11)，可以求得，

(31)

由式(31)求出傅里葉系數(shù).

(33)

(34)

(35)

(36)

根據(jù)位移協(xié)調(diào)條件，位于加筋板板與加強(qiáng)筋分界面位置的質(zhì)點(diǎn)，他們的位移應(yīng)該是相同的.則

ui(x)=w(x,yi),φis(x)=φy(x,yi)

(37)

u0p(x,yi)-φx(x,yi)t/2=u0si(x)+θi(x)h/2

(38)

將式(20)～(22)代入式(13)～(17)，式(33)～(36)代入式(27)～(30)，可以分別用板筋間的耦合力表示廣義位移，聯(lián)合協(xié)調(diào)條件式(37)～式(38)，可以求出耦合力.用矩陣表示為

(39)

(40)

(41)

(42)

式(41)～(42)中的β1=kny1，βN=knyN，η=1-μ.

通過式(39)～(42)就可以求解出板與筋之間的耦合力和力矩，再代入式(20)～(22)，則可以確定加筋板的位移和轉(zhuǎn)角.確定位移之后，求板中應(yīng)力或者筋的應(yīng)力分布，可以根據(jù)幾何方程和本構(gòu)關(guān)系求得，即下面的式(43)～(46).

板的應(yīng)力可以表示為

(43)

(44)

(45)

第i根加強(qiáng)筋的應(yīng)力為

(46)

式中：ω為扇形坐標(biāo).因為加強(qiáng)筋除了彎曲產(chǎn)生的彎曲正應(yīng)力，還有加強(qiáng)筋的約束扭轉(zhuǎn)產(chǎn)生的扇形正應(yīng)力.

3 算例與討論

簡支加筋板長a=2 m、寬b=1.5 m.加強(qiáng)筋為扁鋼，高度h=0.1 m、厚度ts=0.01 m，等間距分布，間距l(xiāng)s=0.5 m，分別布置在y1=0.5 m，y2=1 m的位置.集中力F0=2×104N，作用在板中心的位置，即x0=a/2，y0=b/2.求解加筋板的位移分布和應(yīng)力分布，求解加強(qiáng)筋中的正應(yīng)力分布.

為了驗證本文算法的正確與否，下面采用4種算法，分別命名為方法1～4.方法1采用文獻(xiàn)[10]中的基本原理推導(dǎo)控制微分方程，同時考慮加強(qiáng)筋的扭轉(zhuǎn)變形能，但是求解微分方程的方法不是采用文獻(xiàn)中的最小二乘無網(wǎng)格法，而是采用伽遼金算法.方法2采用ANSYS有限元算法，板采用shell181，共4 800個單元；加強(qiáng)筋采用beam188模擬，共160個單元.方法3 本文算法.方法4 本文方法的退化算法，并忽略整體彎曲中面的偏移構(gòu)件中面的影響.板采用薄板模型，加強(qiáng)筋采用歐拉-伯努利梁模型，則式(1)～(5)退化成一個關(guān)于橫向位移w(x,y)的算式(47).

(47)

加強(qiáng)筋的彎曲方程退化為

(48)

邊界條件和方程的求解也可以仿照前面的方法來調(diào)整.

相應(yīng)板的應(yīng)力計算式，則簡化為

(49)

(50)

加強(qiáng)筋的應(yīng)力計算式修改為

(51)

選取幾個典型部位的位移分布和應(yīng)力分布比較上述4種計算方法的精準(zhǔn)度.圖3為y=b/4，b/2和x=a/4，a/2位置的位移分布；圖4為y=b/4位置處板上表面的x向應(yīng)力分布(以下簡稱x向應(yīng)力分布)；圖5為x=a/4位置處板上表面的y向應(yīng)力分布(以下簡稱y向應(yīng)力分布).圖6為加強(qiáng)筋下邊緣直角點(diǎn)位置處的正應(yīng)力分布(以下簡稱正應(yīng)力分布)，因為外力和加強(qiáng)筋都是關(guān)于y=b/2對稱，所以兩根筋的正應(yīng)力分布完全相同.

圖3 不同位置的位移分布

圖4 y=b/4位置處的x向應(yīng)力分布

圖5 x=a/4位置處的y向應(yīng)力分布

圖6 加強(qiáng)筋中的正應(yīng)力分布

由圖6可見：前3種方法關(guān)于位移的計算結(jié)果都比較接近，僅在中心的位置有細(xì)微差別.板中應(yīng)力分布，方法2和方法3數(shù)據(jù)點(diǎn)比較一致，差別較小，方法1和方法4有的位置與方法2吻合很好，有的位置差別明顯.

有限元法(方法2)是目前普遍認(rèn)同的準(zhǔn)確算法，以方法2的計算值為標(biāo)準(zhǔn)值，則本文方法(方法3)關(guān)于位移的最大計算誤差為3.5%，板的應(yīng)力最大計算誤差為-13.5%，加強(qiáng)筋中應(yīng)力最大計算誤差為-0.1%.以上僅是以圖3～6少數(shù)幾個特殊位置計算結(jié)果得到的誤差分析結(jié)論，但是并非僅這些點(diǎn)的位置吻合較好，其他位置吻合就吻合不好.通過測試比較，其他位置的誤差也在上述范圍之內(nèi).需要說明的是，不能試圖比較載荷作用位置(算例中的板中心點(diǎn))附近的應(yīng)力計算差別，因為在載荷作用位置附近，質(zhì)點(diǎn)的應(yīng)力是無限大，各種方法都不收斂，計算差別很大.這一點(diǎn)，與薄板彎曲的結(jié)論是類似的.

4種計算方法中，方法4數(shù)值處理上最為簡單，精度上在y=b/4位置處的x向應(yīng)力、加強(qiáng)筋中正應(yīng)力與文中方法上差別較大.因為這種方法忽略了彎曲中面偏離加強(qiáng)筋中面的影響，即忽略了板筋之間的剪切作用；還忽略了加強(qiáng)筋的約束扭轉(zhuǎn).加強(qiáng)筋的高度往往比板厚大得多，板筋之間的剪切作用比較明顯，彎曲中面應(yīng)該在靠近板的一側(cè)，忽略這些，將導(dǎo)致加強(qiáng)筋附近板的應(yīng)力和加強(qiáng)筋中應(yīng)力計算誤差較大.

在數(shù)值處理上，本文方法不需要采用諸如伽遼金法或者李茲法求解微分方程，僅需要求解矩陣的逆矩陣或者解線性代數(shù)方程組，編程相對簡單，計算速度明顯比伽遼金算法快.上述4種方法，除方法2是由ANSYS軟件默認(rèn)的迭代次數(shù)和計算精度外，其他3中方法中的三角級數(shù)的迭代累加次數(shù)都不可能是無限次數(shù)的，都是通過編程控制的.雙傅立葉級數(shù)收斂較慢，編程中所有廣義位移表達(dá)式級數(shù)的累加次數(shù)m=n≥11之后則收斂，方法3中Dm矩陣各個元素是單傅立葉級數(shù)，收斂很快，n≥7則收斂.上述計算結(jié)果和誤差分析，都是在級數(shù)收斂這個前提下所完成的.當(dāng)然，累加次數(shù)越多，精度越高，但是計算機(jī)的耗時和內(nèi)存消耗量也是巨大的.

關(guān)于本文方法的適用范圍，需要說明兩點(diǎn).

1)公式推導(dǎo)中，以集中力作為代表，引入狄拉克函數(shù)來表達(dá)的，對于多個集中力的情況，僅需同時列出多個狄拉克函數(shù)累加，可以相應(yīng)推導(dǎo)計算公式；均布載荷也是普遍存在的載荷形式，則需要把推導(dǎo)過程中的Fmn作如下替換即可.

2)上述插圖和算例中是以扁鋼作為代表的，對于其他形式的加強(qiáng)筋，如工字鋼，方鋼，倒T型鋼等也同樣適用，僅需要注意計算加強(qiáng)筋扭轉(zhuǎn)慣性矩應(yīng)該采用下面的公式.

這一點(diǎn)是用能量法所不能達(dá)到的，因為計算應(yīng)變能的時候，不同的加強(qiáng)筋截面需要采用不同的計算方法，如工字鋼、倒T型鋼和方鋼就要分別采用不同的公式計算應(yīng)變能.

4 結(jié) 論

在推導(dǎo)加筋板中的板平衡微分方程，采用了一階剪切變形理論，加強(qiáng)筋(梁)采用了鐵木辛科梁理論，提出了一種“板筋分離法”計算加筋板的橫向彎曲.板與筋之間的耦合力和力矩作為中間變量，加強(qiáng)筋受到板對其的作用力和作用力矩，板也受到加強(qiáng)筋的反作用，以這兩對中間變量代入各自的平衡方程中，求出各自的位移.板與筋共同區(qū)域交界面上，二者具有相同的橫向位移和轉(zhuǎn)角，這就是位移協(xié)調(diào)條件.通過位移協(xié)調(diào)條件求解出中間變量耦合力和力矩，最后求出位移表達(dá)式.

1)通過算例驗證了本文計算方法的正確性，與有限元法相比，位移的計算結(jié)果最大誤差為3.5%；板應(yīng)力計算結(jié)果與有限元法大部分區(qū)域吻合較好，最大誤差-13.5%；加強(qiáng)筋中正應(yīng)力的最大誤差為-0.1%.造成這些局部誤差的主要原因是加強(qiáng)筋的約束扭轉(zhuǎn)沒有完全考慮.

2)文中的計算方法適用性上，可以用于計算集中載荷或均布載荷作用下的加筋板，加強(qiáng)筋的形式可以是扁鋼，工字鋼，方鋼，倒T型鋼等.