薄板彎曲問題有限元法課件

1、 本章將討論彈性板彎曲的有限單元法。當平板的厚度h遠小于其長度a與寬度 時，稱為薄板。對于薄板小撓度問題 ，它的變形完全由橫向撓度w所確定。因此，可以取w和它的若干階導數(shù)作為結(jié)點參數(shù)建立平板單元。目前已經(jīng)提出了非常多的平板單元，但是這里將著重介紹比較常用的矩形單元和一種三角形單元。顯然都不是完全協(xié)調(diào)的單元，但是所得到的計算結(jié)果表明，它們的收斂性和精確度是良好的。 薄板小撓度彎曲問題可視為薄膜彎曲問題，即假設1) Kirchhoff直線法假設。2) 。3) 中面不產(chǎn)生應變。第七章 薄板彎曲問題有限元法第1頁，共81頁。 如圖所示的薄板，取右手坐標系oxyz，使坐標平面oxy位于板的中面，根據(jù)假設

2、2)知：w僅為x、y的函數(shù)，而與z無關，即w = w ( x , y ) 根據(jù)假設1) 得 ，即圖 7-17-1 矩形單元第2頁，共81頁。得上面兩式分別對z積分，并注意 ，即與z無關，得式中 和 是x，y的任意函數(shù)。第3頁，共81頁。 根據(jù)假設3)，可得 ，得而 w = w ( x , y ) (7-32)式中u，v和w是板內(nèi)某點對于坐標軸方向的位移分量。從上面二式可以看出，在平板中面各點u = v = 0，它不產(chǎn)生平面方向的位移，也就是中面不伸長。同時，平板中面的撓度w可以表示板內(nèi)各點的撓度，因為它和坐標z無關。 (7-31)第4頁，共81頁。 利用幾何方程，可以得到板內(nèi)各點的應變分量是(

3、7-33)第5頁，共81頁。 根據(jù)薄板的簡化假定，我們可以把 略去不計，于是板內(nèi)各點的應力可以用撓度表示為式中 (7-35) (7-34)是平板的彈性矩陣，它和平面應力問題中的彈性矩陣完全相同。第6頁，共81頁。 從平板理論知道，若取微元hdxdy，那么在微元上作用著彎矩Mx，My和扭矩Mxy；它是由正應力 和剪應力 在板截面上的合力矩。如果Mx，My和Mxy表示單位寬度上的內(nèi)力矩，于是有式中h是平板厚度。內(nèi)力矩的正方向如圖。(7-36)第7頁，共81頁。 比較(7-34)式和(7-36)式，可以得到用內(nèi)力矩表示的平板應力特別是在平板的上下表面處應力為最大，它是 由以上各式可以看到，平板中面撓

4、度w可以作為基本未知量。如果撓度w為已知，則板中位移、內(nèi)力和應力均可按照上述公式計算。(7-37)第8頁，共81頁。 下面開始講述平板彎曲的有限單元法。 將平板中面用一系列矩形單元劃分，得到一個離散的系統(tǒng)以代替原來的平板，欲使各單元至少在結(jié)點上有撓度及其斜率的連續(xù)性，必須把撓度及其在x和y方向的一階偏導數(shù)指定為結(jié)點位移（或稱廣義位移）。通常將結(jié)點i的位移列陣寫成(7-39)一矩形單元的位移模式zyx第9頁，共81頁。 與之相對應的結(jié)點力列陣可以表示為 它們的符號規(guī)定：對于撓度w和與之對應的結(jié)點力W以沿z軸的正方向為正；對于轉(zhuǎn)角 和與之對應的結(jié)點力矩 ， 則按右手定則標出的矢量沿坐標軸正方向為正

5、。圖7-1中標出的位移和力的方向均為正。 (7-39)第10頁，共81頁。 對于矩形單元，如平面問題中引入一個自然坐標系 來研究單元特性。由于矩形單元的每個結(jié)點有三個位移分量，一個單元有四個結(jié)點共有十二個結(jié)點位移分量，因此我們選取含有十二個參數(shù)的多項式作為位移模式，即(7-40)第11頁，共81頁。 最后兩項的選取是使在單元邊界有三次式的形式。按照上式可以算出轉(zhuǎn)角為 (7-41)第12頁，共81頁。 將矩形單元的四個結(jié)點坐標 分別代入(7-40)式和(7-41)式，就可以得到用十二個參數(shù)表示結(jié)點位移分量的聯(lián)立方程組，求解這十二個方程，從中解出a1至a12再代入(7-40)式，經(jīng)歸并整理后就可以

6、改寫成如下形式或者寫成標準形式其中 (7-43)(7-42)第13頁，共81頁。如果把形函數(shù)寫成通式于是其中記號 和 分別是 ， 。(c) (7-44)第14頁，共81頁。 由(7-31)式可以看到，整個薄板的位移完全由平面在z方向的撓度w所決定，而在中面各點不產(chǎn)生x和y方向位移。因此薄板所可能產(chǎn)生的剛性位移就只有沿z方向的平動以及繞x和y軸的轉(zhuǎn)動，而對于z軸方向的旋轉(zhuǎn)是沒有的。位移模式(7-40)式中是前三項反映了薄板單元的這三個剛體位移。再由(7-33)式看到，板內(nèi)各點的應變完全由撓度w的三個二階導數(shù)所決定。如果這三個二階導數(shù)不隨坐標而變化，則描述平板單元的一個常應變狀態(tài)，(7-40)式中

7、的第四、五、六三個二次項反映了這個常應變狀態(tài)（或稱常曲率狀態(tài)）。因此，我們總是能夠保證存在一組結(jié)點位移，可以反映單元的剛體位移和常應變狀態(tài)，因此，這個矩形單元是完備的。 從(7-33)式和(7-34)式可以看出，應變和應力是有撓度w的二階偏導數(shù)所決定。因此，如果要得到一個協(xié)調(diào)的單元還要求在單元的交界面上有斜率的連續(xù)性，這個要求經(jīng)常使問題復雜化。第15頁，共81頁。 由(7-42)和(7-44)式可以看出，在單元邊界上撓度和撓度沿切線邊界方向的偏導數(shù)，可以通過邊界上的結(jié)點位移所唯一地決定，但是撓度沿邊界法線方向的偏導數(shù)則不然，也就是說，w和 的值在單元交界線之間是連續(xù)的，而對于 卻不連續(xù)；s表示

8、交界線切線方向而n表示交界線法線方向，因此，我們現(xiàn)在所討論的單元是非協(xié)調(diào)元，或稱為不完全協(xié)調(diào)單元。 以 的ij邊界為例說明ssn1n2ij第16頁，共81頁。 該邊界上兩端點i , j共有4個已知條件：將這4個條件代入 中，就可以完全確定4個常數(shù)c1,c2,c3,c4。如果 邊界 是兩相鄰單元的公共邊界，則兩個單元分別按上述4 個條件所確定的常數(shù)c1,c2,c3,c4也一定相同，即兩相鄰單元的公共邊界、上有相同的撓度w。這表明，所選取的位移模式w滿足了相鄰單元的撓度在公共邊界上的連續(xù)條件。第17頁，共81頁。再由式(7-41)的第一式看出，在單元 邊界上的法線轉(zhuǎn)角 也是x（或 ）的三次多項式上

9、式仍需要兩端點i , j有4個已知條件來確定常數(shù)d1,d2,d3,d4，但是，現(xiàn)在只有 和 兩個條件，不可能確定出4個常數(shù)d1,d2,d3,d4。因此，板單元整個 公共邊界上的法線轉(zhuǎn)角 是不連續(xù)的，只有在公共邊界的兩端點i , j上有共同的法線轉(zhuǎn)角。第18頁，共81頁。 將(7-42)式代入幾何方程式(7-33)，可以將單元應變用結(jié)點位移列陣表示為式中記號 等分別表示 。 (7-45)(7-46)二矩形單元的剛度矩陣第19頁，共81頁。 按照(c)式和(7-44)式可以算出(d)第20頁，共81頁。 于是單元剛度矩陣可以寫成如下形式其中子矩陣的計算公式是(7-47)第21頁，共81頁。 把(7

10、-35)式和(7-36)式代入上式，并完成對z的積分，于是有式中 它就是彈性薄板的彎曲剛度。如果再利用(d)式把(7-48)式展開并完成全部積分，就可以得到子矩陣(7-48) (7-49) (7-50)第22頁，共81頁。式中的九個元素的顯式如下第23頁，共81頁。(7-51) 式中第24頁，共81頁。 如果平板單元受有分布橫向載荷q的作用，于是等效結(jié)點力是 當q = q0為常量時，將(7-14)式代入上式并進行積分，于是得 (7-52)三矩形單元的等效結(jié)點力和內(nèi)力矩計算第25頁，共81頁。 最后，由(7-37)式知道，若要計算平板應力列陣 ，必需算出內(nèi)力矩列陣 。而對于 的計算，只要在(7-

11、33)式和(7-45)式中求得第26頁，共81頁。 再把上式代入(7-36)式中，可以得到式中 (7-53)第27頁，共81頁。 若將(d)式代入上式，則得 第28頁，共81頁。 由于矩形單元在使用上受到平板形狀的限制，而采用三角形單元可以較好地反映邊界形狀。 根據(jù)板單元每個結(jié)點三個位移，而三角形單元三個結(jié)點，于是被采用的位移模式應該包含9個參數(shù)，而x和y的完全三次多項式共計十項。若以它為基礎構(gòu)造位移模式，必須在其中刪去一項。而三次方項刪去任何一項，都不能保持對于x和y的對稱性，有人建議取7-2 三角形單元一三角形單元的位移模式第29頁，共81頁。以達到減少一個待定系數(shù)并保持對稱性的目的?？上?/p>

12、在此情況下，對于二個邊界分別平行于x軸和y軸的等腰三角形單元，確定的代數(shù)方程系數(shù)矩 是奇異的，因此陣不能確定，此方案不行。還有另一種方案是將單元中心撓度w也作為一個參數(shù)，但按此方案導出的單元是不收斂的。因此，在直角坐標系中構(gòu)造三角形板單元的撓度插值函數(shù)是困難的，而在面積坐標下進行這項工作可行。 用1、2、3代替i、j、m，則面積坐標的一次、二次、三次式分別有以下各項。 一次 二次 三次 第30頁，共81頁。 容易看出三次式的最后一項 （注意L3 = 1L1L2）本身和它的兩個一階偏導數(shù)，在三個角點處的值等于零，對于確定9個參數(shù)無用，因此自然可以刪去而利用前面九項來構(gòu)造位移模式。但是，由這個不完

13、全的三次多項式構(gòu)成的位移模式，不能保證有獨立的線性項和二次項；也就是說，剛體位移和常應變準則，可能不被滿足。為了這一點，可假設位移模式是 (7-25)第31頁，共81頁。式中前三項反映剛體位移，次三項對應于常應變。二次項只取了后三項是為了用結(jié)點位移表示參數(shù) 時考慮計算上的方便。同理，三次項不取前三項，剩下六個，只能挑選三個或進行某種線性組合。為了考慮每項面積坐標對稱地出現(xiàn)，作出了如上的最簡單可行的線性組合。組合未取“+”號，是由于，所以使用加號的最簡單線性組合是不合宜的。 為了將位移模式寫成標準形式，就需要求得形函數(shù)。為了方便起見，求形函數(shù)的工作可以分成兩步進行。第一步是選取w、 、 作為結(jié)點

14、自由度，求得對于它們的形函數(shù)，在這里把L3 = 1L1L2看作是L1和L2的函數(shù)。第二步是利用關系式第32頁，共81頁。 將第一步中所用的結(jié)點自由度變換成(7-38)式所指定的結(jié)點位移，再通過合并整理就很容易地得到形函數(shù)。式中b1 = y2 y3，c1 = x3 x1。對于b2、c2的值可以用下標輪換定出。 (7-26)第33頁，共81頁。 現(xiàn)在來決定參數(shù) 到 。將三角形單元的三個結(jié)點的面積坐標代入(7-25)式，立即得到 ， ， 。利用(7-25)式計算 和 ，得到第34頁，共81頁。 將結(jié)點的面積坐標代入上式，得六個方程如下第35頁，共81頁。式中表示對Li的偏導數(shù)在j點的值。從上式解得第

15、36頁，共81頁。將上式代入(7-25)式，并歸并 和 前的各項，就可以得到對應于它們的形函數(shù)。 利用(7-26)式，將 和 變換為 和 ，于是得到相應于 和 的形函數(shù) 和 的計算公式。第37頁，共81頁。 最后得如下形式的形函數(shù)為(d)第38頁，共81頁。用下標輪換可得結(jié)點2和3的形函數(shù)，將上式寫成矩陣形式其中 (7-29)而 是一個的系數(shù)矩陣，它是第39頁，共81頁。式中它們分別是結(jié)點1、2和3的面積坐標。位移模式可寫成如下的標準形式第40頁，共81頁。 可以驗證在相鄰單元間的撓度是連續(xù)的，但它的法向斜率仍不連續(xù)。事實上，在任何一條邊上，撓度可表達成邊線方向s的三次式，并且不包含與此邊相對

16、結(jié)點的結(jié)點位移（因為與它相對應的形函數(shù)在此邊上等于零）。也就是說，一條邊上的撓度可以由端部兩個結(jié)點處的w和 所完全決定，而對于 則不然。因此這個三角形單元是一個完備的非協(xié)調(diào)單元。第41頁，共81頁。 在推導剛度矩陣和彎矩公式時，要計算 、 和 由于形函數(shù)是用面積坐標表示的，因此必須寫出兩個坐標系中的偏導數(shù)之間的關系。仍然取L1、L2作為獨立坐標，而L3 = 1L1L2作為L1和L2的函數(shù)，根據(jù)復合函數(shù)的求導數(shù)規(guī)則，并且利用坐標變換公式，可以得到下列兩個關系式。二三角形單元的剛度矩陣第42頁，共81頁。式中 是三角形的面積，而 將w的標準式代入幾何方程，得單元應變列陣第43頁，共81頁。式中記號

17、 N i,11等表示 N i對于L1的兩次偏導數(shù)等。將(7-29)式代入上式，得到第44頁，共81頁。由于 L 是面積坐標的三次函數(shù)，它的三個對L1和L2的二階偏導數(shù)將是面積坐標的一次函數(shù)。把這三個偏導數(shù)算出后，容易把上式寫成如下形式其中第45頁，共81頁。(7-38)第46頁，共81頁。而可以將 C Ai 乘出并記作矩陣 Gi ，于是式中(7-40)第47頁，共81頁。 若將單元剛度矩陣寫成如下形式第48頁，共81頁。其中子矩陣若令(7-41)第49頁，共81頁。于是得式中D是平板的彎曲剛度。 注意到(7-38)式，上式的積分是容易計算的。實際上(7-43) (7-42)第50頁，共81頁。

18、式中它的積分是而第51頁，共81頁。 把(7-40)和(7-43)式代入(7-42)式并把它展開，可得單元剛度矩陣的子矩陣第52頁，共81頁。式中 H 的所有元素，按(7-41)式展開進行計算。第53頁，共81頁。 如果平板單元受有分布橫向載荷q的作用，于是等效結(jié)點力是這里以使用(d)式所示的形函數(shù)較為方便。三三角形單元的等效結(jié)點力和內(nèi)力矩的計算第54頁，共81頁。 當q = q0是常量時，將(d)式代入上式，并積分得第55頁，共81頁。 對于內(nèi)力矩陣，按下式計算式中(7-49)第56頁，共81頁。把上式代入(7-49)式并進行一系列運算之后，可以得到內(nèi)力矩列陣式中 對于任意形狀平板，它的邊界

19、條件可能是指定為沿著曲線邊界切線方向的彎矩Ms或轉(zhuǎn)角 和沿著法線方向的扭矩Mn或轉(zhuǎn)角 ，這里所用的記號與一般平板理論書中的記號恰好相反。這里取邊界的外法線方向n為正方向，而使用右手坐標系定出切線s的正方向。于是，有下列關系式存在第57頁，共81頁。式中c、s是外法線n方向?qū)τ趚和y軸的方向余弦。 利用上式可以對剛度矩陣作某種變動，使之達到結(jié)點位移和結(jié)點力的變換，從而直接利用曲線邊界上的邊界條件。這種變動可以在單元剛度矩陣中進行，也可以在整體剛度矩陣中進行。第58頁，共81頁。 上述矩形單元和三角形薄板單元都是非協(xié)調(diào)元，且都會產(chǎn)生在單元公共邊上的法線轉(zhuǎn)角 （或法向斜率 ）不連續(xù)問題。為了實現(xiàn)薄板

20、單元的協(xié)調(diào)性（協(xié)調(diào)板元或保續(xù)板元），完善板元的計算理論和提高計算精度，人們做了大量研究，提出了很多方法。其中一種方法是增加單元結(jié)點自由度數(shù)目，例如在矩形板元的每個結(jié)點上增加一個扭率 ，使單元變成16個自由度的協(xié)調(diào)板單元；在三角形板單元中每個結(jié)點上增加 、 、 做為結(jié)點自由度，且在每邊中點取其法線斜率 做為結(jié)點自由度、而構(gòu)造出21個自由度的協(xié)調(diào)三角形板單元。 7-3 考慮橫向剪切影響的平板彎曲單元第59頁，共81頁。 另一種方法，是將三角形（或四邊形）板單元劃分為3個子三角形，把每個三角形各邊中點的法線斜率 做為自由度，每個子三角形都有12個自由度，使各子三角形板元之間協(xié)調(diào)，然后再根據(jù)原三角形板

21、元內(nèi)部的連續(xù)性和限制條件，利用“凝聚法”消去內(nèi)自由度，從而構(gòu)造出協(xié)調(diào)的三角形板元 （12個外自由度，3個內(nèi)自由度）或協(xié)調(diào)的四邊形板元。第60頁，共81頁。 第一種方法的明顯缺點：在實際應用時涉及到高階導數(shù)的自由度的邊界條件，難以處理；它不是一種普遍適用的方法，在某些板問題中，曲率或扭率在結(jié)點上不一定連續(xù)，在有限元計算中強令其連續(xù)，當然會使所得結(jié)果不可能收斂到精確解。 第二種方法在SAP-5程序中得到應用。但該法的計算公式比較復雜，消去每個單元的內(nèi)自由度的凝聚過程所耗計算時間比較長。 在板彎曲問題的有限元法中，構(gòu)造協(xié)調(diào)元的困難是單元公共邊上的法向轉(zhuǎn)角 的連續(xù)性難以滿足，如果考慮板橫向剪切變形的影

22、響，放棄經(jīng)典薄板理論中的中面法線n-n始終保持為直線的假設，就可以繞開這個困難，而使板問題的有限元分析前進一步。第61頁，共81頁。 一般來說，板變形前的中面法線n-n，在變形后將變成一條曲線，但可以近似地用一條直線n2-n2表示，n2-n2線已不是經(jīng)典薄板理論中變形后的中面法線n-n。此時，n2-n2線繞x軸和y軸的轉(zhuǎn)角仍用 和 表示，但 。 經(jīng)典薄板理論中的其它兩個假設，仍然有效。第62頁，共81頁。 根據(jù)上述假定，板內(nèi)任意點的3 個位移分量具有如下形式把上式代入幾何方程(6-33)式，可以得到平板應變列陣是(7-54)(7-53)第63頁，共81頁。 通過應力應變關系，可以得到應力列陣如

23、下式中彈性矩陣是而它的子矩陣(a) (7-55) (7-56)第64頁，共81頁。 由(7-53)式看到，平板的變形完全是由中面撓度w及其法線繞x軸和y軸兩個轉(zhuǎn)角 和 所決定。在每個結(jié)點上取它們作為自由度，構(gòu)造一個八結(jié)點平板單元，它與第三章的八結(jié)點等參元密切相關。 對于八結(jié)點平板單元，要確定它的形狀除了單元中面的形狀外，還需要知道單元厚度h。中面的形狀可以通過坐標變換式(3-2)由八個結(jié)點坐標來確定，而單元的厚度可以利用形函數(shù)表達式(3-1)通過八個結(jié)點處的厚度近似地用插入法確定出。顯然，自然坐標系是位于板中面內(nèi)的一種曲線坐標。第65頁，共81頁。 中面上任意點的撓度w以及法線轉(zhuǎn)角 和 ，同樣

24、可以利用(3-1)式所表示的形函數(shù)由結(jié)點值進行插值得到，即代入(7-53)式，于是得到位移模式如下容易看出單元是協(xié)調(diào)的。這里沒有必要引進第三個自然坐標 ，因為引進它并沒有帶來任何方便，顯然 = 2z/h。 (b)第66頁，共81頁。 將(b)式代入(7-54)式?？梢缘玫接媒Y(jié)點位移列陣表示的應變分量式中(7-59)(7-58)第67頁，共81頁。而 (7-60)第68頁，共81頁。 將(7-58)式代入(7-55)式，可以作出用結(jié)點位移列陣表示的應力分量式中而第69頁，共81頁。由此按定義可以計算內(nèi)力，它的表達式如下第70頁，共81頁。式中 (7-63) 將單元剛度矩陣寫成如(3-12)式的形式，其中子矩陣可以按下式計算(7-64)第71頁，共81頁。若命于是，(5-64)式可以寫成第72頁，共81頁。式中 是雅可比行