彈性力學(xué)平面問題的有限元法(共30頁(yè))

1、精選優(yōu)質(zhì)文檔-傾情為你奉上Mmm3 彈性力學(xué)平面問題的有限元法本章包括以下的內(nèi)容：3.1彈性力學(xué)平面問題的基本方程3.2單元位移函數(shù)3.3單元載荷移置3.4單元?jiǎng)偠染仃?.5單元?jiǎng)偠染仃嚨男再|(zhì)與物理意義3.6整體分析3.7約束條件的處理3.8整體剛度矩陣的特點(diǎn)與存儲(chǔ)方法3.9方程組解法3.1彈性力學(xué)平面問題的基本方程彈性力學(xué)是研究彈性體在約束和外載荷作用下應(yīng)力和變形分布規(guī)律的一門學(xué)科。在彈性力學(xué)中針對(duì)微小的單元體建立基本方程，把復(fù)雜形狀彈性體的受力和變形分析問題歸結(jié)為偏微分方程組的邊值問題。彈性力學(xué)的基本方程包括平衡方程、幾何方程、物理方程。彈性力學(xué)的基本假定如下：1）完全彈性，2）連續(xù)，3）

2、均勻，4）各向同性，5）小變形。3.1.1基本變量彈性力學(xué)中的基本變量為體力、面力、應(yīng)力、位移、應(yīng)變，各自的定義如下。體力體力是分布在物體體積內(nèi)的力，例如重力和慣性力。面力面力是分布在物體表面上的力，例如接觸壓力、流體壓力。應(yīng)力 物體受到約束和外力作用，其內(nèi)部將產(chǎn)生內(nèi)力。物體內(nèi)某一點(diǎn)的內(nèi)力就是應(yīng)力。圖3.1如圖3.1假想用通過物體內(nèi)任意一點(diǎn)p的一個(gè)截面mn將物理分為、兩部分。將部分撇開，根據(jù)力的平衡原則，部分將在截面mn上作用一定的內(nèi)力。在mn截面上取包含p點(diǎn)的微小面積，作用于面積上的內(nèi)力為。令無限減小而趨于p點(diǎn)時(shí)，的極限S就是物體在p點(diǎn)的應(yīng)力。 應(yīng)力S在其作用截面上的法向分量稱為正應(yīng)力，用表

3、示；在作用截面上的切向分量稱為剪應(yīng)力，用表示。顯然，點(diǎn)p在不同截面上的應(yīng)力是不同的。為分析點(diǎn)p的應(yīng)力狀態(tài)，即通過p點(diǎn)的各個(gè)截面上的應(yīng)力的大小和方向，在p點(diǎn)取出的一個(gè)平行六面體，六面體的各楞邊平行于坐標(biāo)軸。圖3.2將每個(gè)上的應(yīng)力分解為一個(gè)正應(yīng)力和兩個(gè)剪應(yīng)力，分別與三個(gè)坐標(biāo)軸平行。用六面體表面的應(yīng)力分量來表示p點(diǎn)的應(yīng)力狀態(tài)。應(yīng)力分量的下標(biāo)約定如下：第一個(gè)下標(biāo)表示應(yīng)力的作用面，第二個(gè)下標(biāo)表示應(yīng)力的作用方向。，第一個(gè)下標(biāo)x表示剪應(yīng)力作用在垂直于X軸的面上，第二個(gè)下標(biāo)y表示剪應(yīng)力指向Y軸方向。正應(yīng)力由于作用表面與作用方向垂直，用一個(gè)下標(biāo)。表示正應(yīng)力作用于垂直于X軸的面上，指向X軸方向。 應(yīng)力分量的方向定

4、義如下： 如果某截面上的外法線是沿坐標(biāo)軸的正方向，這個(gè)截面上的應(yīng)力分量以沿坐標(biāo)軸正方向?yàn)檎?；如果某截面上的外法線是沿坐標(biāo)軸的負(fù)方向，這個(gè)截面上的應(yīng)力分量以沿坐標(biāo)軸負(fù)方向?yàn)檎?。剪?yīng)力互等：物體內(nèi)任意一點(diǎn)的應(yīng)力狀態(tài)可以用六個(gè)獨(dú)立的應(yīng)力分量、來表示。位移位移就是位置的移動(dòng)。物體內(nèi)任意一點(diǎn)的位移，用位移在x，y，z坐標(biāo)軸上的投影u、v、w表示。應(yīng)變物體的形狀改變可以歸結(jié)為長(zhǎng)度和角度的改變。各線段的單位長(zhǎng)度的伸縮，稱為正應(yīng)變，用表示。兩個(gè)垂直線段之間的直角的改變，用弧度表示，稱為剪應(yīng)變，用表示。物體內(nèi)任意一點(diǎn)的變形，可以用六個(gè)應(yīng)變分量表示。3.1.2平面應(yīng)力和平面應(yīng)變問題彈性體在滿足一定條件時(shí)，其變形和

5、應(yīng)力的分布規(guī)律可以用在某一平面內(nèi)的變形和應(yīng)力的分布規(guī)律來代替，這類問題稱為平面問題。平面問題分為平面應(yīng)力問題和平面應(yīng)變問題。1） 平面應(yīng)力問題設(shè)有很薄的等厚薄板，只在板邊上受到平行于板面并且不沿厚度變化的面力，體力也平行于板面且不沿厚度變化。圖3.3設(shè)板的厚度為t，在板面上：， ， 由于平板很薄，外力不沿厚度變化，因此在整塊板上有， ， 剩下平行于XY平面的三個(gè)應(yīng)力分量未知。2）平面應(yīng)變問題設(shè)有很長(zhǎng)的柱形體，支承情況不沿長(zhǎng)度變化，在柱面上受到平行于橫截面而且不沿長(zhǎng)度變化的面力，體力也如此分布。圖3.4以柱體的任一橫截面為XY平面，任一縱線為Z軸。假定該柱體為無限長(zhǎng)，則任一截面都可以看作對(duì)稱面。

6、由對(duì)稱性，由于沒有Z方向的位移，Z方向的應(yīng)變。未知量為平行于XY平面的三個(gè)應(yīng)力分量，物體在Z方向處于自平衡狀態(tài)。3.1.3平衡方程彈性力學(xué)中，在物體中取出一個(gè)微小單元體建立平衡方程。平衡方程代表了力的平衡關(guān)系，建立了應(yīng)力分量和體力分量之間的關(guān)系。對(duì)于平面問題，在物體內(nèi)的任意一點(diǎn)有，（3-1）3.1.4幾何方程由幾何方程可以得到位移和變形之間的關(guān)系。對(duì)于平面問題，在物體內(nèi)的任意一點(diǎn)有，（3-2）剛體位移由位移u=0，v=0可以得到應(yīng)變分量為零，反過來，應(yīng)變分量為零則位移分量不為零。應(yīng)變分量為零時(shí)的位移稱為剛體位移。剛體位移代表了物體在平面內(nèi)的移動(dòng)和轉(zhuǎn)動(dòng)。由可以得到剛體位移為以下形式，由可得，將代

7、入可得，積分后得到，得到位移分量，當(dāng)時(shí)，物體內(nèi)任意一點(diǎn)都沿x方向移動(dòng)相同的距離，可見代表物體在x方向上的剛體平移。當(dāng)時(shí)，物體內(nèi)任意一點(diǎn)都沿y方向移動(dòng)相同的距離，可見代表物體在y方向上的剛體平移。當(dāng)時(shí)，可以假定，此時(shí)的物體內(nèi)任意一點(diǎn)P（x，y）的位移分量為，P點(diǎn)位移與y軸的夾角為，P點(diǎn)合成位移為，r為P點(diǎn)到原點(diǎn)的距離，可見代表物體繞z軸的剛體轉(zhuǎn)動(dòng)。3.1.5物理方程彈性力學(xué)平面問題的物理方程由廣義虎克定律得到。1）平面應(yīng)力問題的物理方程（3-3）平面應(yīng)力問題有，2）平面應(yīng)變問題的物理方程（3-4）平面應(yīng)變問題有，在平面應(yīng)力問題的物理方程中，將E替換為、替換為，可以得到平面應(yīng)變問題的物理方程；在平

8、面應(yīng)變問題的物理方程中，將E替換為、替換為，可以得到平面應(yīng)力問題的物理方程。圖3.5求解彈性力學(xué)平面問題，可以歸結(jié)為在任意形狀的平面區(qū)域內(nèi)已知控制方程、在位移邊界上約束已知、在應(yīng)力邊界上受力條件已知的邊值問題。然后以應(yīng)力分量為基本未知量求解，或以位移作為基本未知量求解。如果以位移作為未知量求解，求出位移后，由幾何方程可以計(jì)算出應(yīng)變分量，得到物體的變形情況；再由物理方程計(jì)算出應(yīng)力分量，得到物體的內(nèi)力分布，就完成了對(duì)彈性力學(xué)平面問題的分析。3.2單元位移函數(shù)根據(jù)有限元法的基本思路，將彈性體離散成有限個(gè)單元體的組合，以結(jié)點(diǎn)的位移作為未知量。彈性體內(nèi)實(shí)際的位移分布可以用單元內(nèi)的位移分布函數(shù)來分塊近似地

9、表示。在單元內(nèi)的位移變化可以假定一個(gè)函數(shù)來表示，這個(gè)函數(shù)稱為單元位移函數(shù)、或單元位移模式。對(duì)于彈性力學(xué)平面問題，單元位移函數(shù)可以用多項(xiàng)式表示，（3-5）多項(xiàng)式中包含的項(xiàng)數(shù)越多，就越接近實(shí)際的位移分布，越精確。具體取多項(xiàng)，由單元形式來確定。即以結(jié)點(diǎn)位移來確定位移函數(shù)中的待定系數(shù)。圖3.6如圖3.6所示的3結(jié)點(diǎn)三角形單元，結(jié)點(diǎn)I、J、M的坐標(biāo)分別為、，結(jié)點(diǎn)位移分別為、。六個(gè)節(jié)點(diǎn)位移只能確定六個(gè)多項(xiàng)式的系數(shù)，所以3結(jié)點(diǎn)三角形單元的位移函數(shù)如下，（3-6）將3個(gè)結(jié)點(diǎn)上的坐標(biāo)和位移分量代入公式（3-6）就可以將六個(gè)待定系數(shù)用結(jié)點(diǎn)坐標(biāo)和位移分量表示出來。將水平位移分量和結(jié)點(diǎn)坐標(biāo)代入（3-6）中的第一式，寫

10、成矩陣形式，（3-7）令，則有（3-8），A為三角形單元的面積。T的伴隨矩陣為，（3-9）令（3-10）則（3-11）同樣，將垂直位移分量與結(jié)點(diǎn)坐標(biāo)代入公式（3-6）中的第二式，可得，（3-12）將（3-11）、（3-12）代回（3-6）整理后可得，令（下標(biāo)i，j，m輪換）可得（3-13）單元內(nèi)的位移記為單元的結(jié)點(diǎn)位移記為單元內(nèi)的位移函數(shù)可以簡(jiǎn)寫成，（3-14）把N稱為形態(tài)矩陣，Ni稱為形態(tài)函數(shù)。選擇單元位移函數(shù)應(yīng)滿足以下條件：1） 反映單元的剛體位移與常量應(yīng)變。2）相鄰單元在公共邊界上的位移連續(xù)，即單元之間不能重疊，也不能脫離。由（3-6）可以將單元位移表示成以下的形式，反映了剛體位移和常應(yīng)

11、變。 單元位移函數(shù)是線性插值函數(shù)，因此單元邊界上各點(diǎn)的位移可以由兩個(gè)結(jié)點(diǎn)的位移完全確定。兩個(gè)單元的邊界共用兩個(gè)結(jié)點(diǎn)，所以邊界上的位移連續(xù)。 形態(tài)函數(shù)Ni具有以下性質(zhì)：1） 在單元結(jié)點(diǎn)上形態(tài)函數(shù)的值為1或?yàn)?。2）在單元中的任意一點(diǎn)上，三個(gè)形態(tài)函數(shù)之和等于1。用來計(jì)算三角形面積時(shí)，要注意單元結(jié)點(diǎn)的排列順序，當(dāng)三個(gè)結(jié)點(diǎn)i，j，m取逆時(shí)針順序時(shí)，；當(dāng)三個(gè)結(jié)點(diǎn)i，j，m取順時(shí)針順序時(shí)，。例題：如圖3.7所示等腰三角形單元，求其形態(tài)矩陣N。解：由在公式中輪換下標(biāo)可以計(jì)算得，三角形積為形態(tài)函數(shù)為形態(tài)矩陣為三角形面積的計(jì)算公式可得，如果把三個(gè)結(jié)點(diǎn)按順時(shí)針方向排列，即i（a，0），j（0，0），m（0，a）3

12、.3單元載荷移置有限元法的求解對(duì)象是單元的組合體，因此作用在彈性體上的外力，需要移置到相應(yīng)的結(jié)點(diǎn)上成為結(jié)點(diǎn)載荷。載荷移置要滿足靜力等效原則。靜力等效是指原載荷與結(jié)點(diǎn)載荷在任意虛位移上做的虛功相等。單元的虛位移可以用結(jié)點(diǎn)的虛位移表示為，（3-15）令結(jié)點(diǎn)載荷為1）集中力的移置如圖3.7所示，在單元內(nèi)任意一點(diǎn)作用集中力圖3.8由虛功相等可得，由于虛位移是任意的，則（3-16）例題1：在均質(zhì)、等厚的三角形單元ijm的任意一點(diǎn)p（xp，yp）上作用有集中載荷。2） 體力的移置令單元所受的均勻分布體力為 由虛功相等可得，（3-17）3） 分布面力的移置設(shè)在單元的邊上分布有面力，同樣可以得到結(jié)點(diǎn)載荷，（3

13、-18）例題2：設(shè)有均質(zhì)、等厚的三角形單元ijm，受到沿y方向的重力載荷qy的作用。求均布體力移置到各結(jié)點(diǎn)的載荷。同理，例題3：在均質(zhì)、等厚的三角形單元ijm的ij邊上作用有沿x方向按三角形分布的載荷，求移置后的結(jié)點(diǎn)載荷。取局部坐標(biāo)s，在i點(diǎn)s=0，在j點(diǎn)s=l，L為ij邊的長(zhǎng)度。在ij邊上，以局部坐標(biāo)表示的插值函數(shù)為，載荷為3.4單元?jiǎng)偠染仃?根據(jù)單元的位移函數(shù)，由幾何方程可以得到單元的應(yīng)變表達(dá)式，（3-19）記為，B矩陣稱為幾何矩陣。B矩陣可以表示為分塊矩陣的形式（3-20） 由物理方程，可以得到單元的應(yīng)力表達(dá)式，（3-21） D稱為彈性矩陣，對(duì)于平面應(yīng)力問題， 定義為應(yīng)力矩陣。將應(yīng)力矩陣

14、分塊表示為，（3-22）應(yīng)用虛功原理可以建立單元結(jié)點(diǎn)位移與結(jié)點(diǎn)力的關(guān)系矩陣，單元?jiǎng)偠染仃?。虛功原理：在外力作用下處于平衡狀態(tài)的彈性體，如果發(fā)生了虛位移，則所有外力在虛位移上做的虛功等于內(nèi)應(yīng)力在虛應(yīng)變上做的虛功。單元的結(jié)點(diǎn)力記為單元的虛應(yīng)變?yōu)閱卧耐饬μ摴?，單元的?nèi)力虛功為，由虛功原理可得，（3-23）（3-24）定義為單元?jiǎng)偠染仃?。?結(jié)點(diǎn)等厚三角形單元中B和D的分量均為常量，則單元?jiǎng)偠染仃嚳梢员硎緸?，?-25）單元?jiǎng)偠染仃嚤硎緸榉謮K矩陣：（3-26）（3-27）3.5單元?jiǎng)偠染仃嚨男再|(zhì)與物理意義（一）單元?jiǎng)偠染仃嚨奈锢硪饬x假設(shè)單元的結(jié)點(diǎn)位移如下：由，得到結(jié)點(diǎn)力如下：（3-28）表示i結(jié)點(diǎn)

15、在水平方向產(chǎn)生單位位移時(shí)，在結(jié)點(diǎn)i的水平方向上需要施加的結(jié)點(diǎn)力。表示i結(jié)點(diǎn)在水平方向產(chǎn)生單位位移時(shí)，在結(jié)點(diǎn)i的垂直方向上需要施加的結(jié)點(diǎn)力。選擇不同的單元結(jié)點(diǎn)位移，可以得到單元?jiǎng)偠染仃囍忻總€(gè)元素的物理含義：表示s結(jié)點(diǎn)在水平方向產(chǎn)生單位位移時(shí)，在結(jié)點(diǎn)r的水平方向上需要施加的結(jié)點(diǎn)力。表示s結(jié)點(diǎn)在水平方向產(chǎn)生單位位移時(shí)，在結(jié)點(diǎn)r的垂直方向上需要施加的結(jié)點(diǎn)力。表示s結(jié)點(diǎn)在垂直方向產(chǎn)生單位位移時(shí)，在結(jié)點(diǎn)r的水平方向上需要施加的結(jié)點(diǎn)力。表示s結(jié)點(diǎn)在垂直方向產(chǎn)生單位位移時(shí)，在結(jié)點(diǎn)r的垂直方向上需要施加的結(jié)點(diǎn)力。因此單元?jiǎng)偠染仃囍忻總€(gè)元素都可以理解為剛度系數(shù)，即在結(jié)點(diǎn)產(chǎn)生單位位移時(shí)需要施加的力。（二）單元?jiǎng)偠染?/p>

16、陣的性質(zhì)1） 對(duì)稱性利用分塊矩陣的性質(zhì)證明如下：即2） 奇異性即單元?jiǎng)偠染仃嚨男辛惺綖榱?，。將定單元產(chǎn)生了x方向的剛體移動(dòng)，此時(shí)對(duì)應(yīng)的單元結(jié)點(diǎn)力為零?？梢缘玫?，在單元?jiǎng)偠染仃囍?，3，5列中對(duì)應(yīng)行的系數(shù)相加為零，由行列式的性質(zhì)可知，。同樣如果假定單元產(chǎn)生了y方向上的剛體位移，可以得到，在單元?jiǎng)偠染仃囍?，4，6列中對(duì)應(yīng)行的系數(shù)相加為零。3.6整體分析 得到了單元?jiǎng)偠染仃嚭?，要將單元組成一個(gè)整體結(jié)構(gòu)，根據(jù)結(jié)點(diǎn)載荷平衡的原則進(jìn)行分析，即整體分析。整體分析包括以下4個(gè)步驟：1） 建立整體剛度矩陣，2） 根據(jù)支承條件修改整體剛度矩陣，3） 解方程組，求出結(jié)點(diǎn)的位移，4） 根據(jù)結(jié)點(diǎn)位移，求出單元的應(yīng)變和

17、應(yīng)力。在這里把結(jié)點(diǎn)位移作為基本未知量求解。如何得到整體剛度矩陣？基本方法是剛度集成法，即整體剛度矩陣是單元?jiǎng)偠染仃嚨募?。圖3.9 如圖3.9所示，一個(gè)劃分為6個(gè)結(jié)點(diǎn)、4個(gè)單元的結(jié)構(gòu)。得到了每個(gè)單元的單元?jiǎng)偠染仃嚭?，要集成為整體剛度矩陣。3.6.1剛度集成法的物理意義 由單元?jiǎng)偠染仃嚨奈锢硪饬x可知，單元?jiǎng)偠染仃嚨南禂?shù)是由單元結(jié)點(diǎn)產(chǎn)生單位位移時(shí)引起的單元結(jié)點(diǎn)力。在如圖3.9所示的結(jié)構(gòu)中，使結(jié)點(diǎn)3產(chǎn)生單位位移時(shí)，在單元（1）中的結(jié)點(diǎn)2上引起結(jié)點(diǎn)力。由于結(jié)點(diǎn)2、3同時(shí)屬于單元（1）、（3），在單元（2）中的結(jié)點(diǎn)2上同樣也引起結(jié)點(diǎn)力，因此，在整體結(jié)構(gòu)中當(dāng)結(jié)點(diǎn)3產(chǎn)生位移時(shí)，結(jié)點(diǎn)2上的結(jié)點(diǎn)力應(yīng)該是單元（1

18、）、（2）在結(jié)點(diǎn)2上的結(jié)點(diǎn)力的疊加。剛體集成法即結(jié)構(gòu)中的結(jié)點(diǎn)力是相關(guān)單元結(jié)點(diǎn)力的疊加，整體剛度矩陣的系數(shù)是相關(guān)單元的單元?jiǎng)偠染仃囅禂?shù)的集成。結(jié)點(diǎn)3在整體剛度矩陣的對(duì)應(yīng)系數(shù)，應(yīng)該是單元（1）、（3）、（4）中對(duì)應(yīng)系數(shù)的集成。3.6.2剛度矩陣集成的規(guī)則1） 將單元?jiǎng)偠染仃囍械拿總€(gè)分塊放到在整體剛度矩陣中的對(duì)應(yīng)位置上，得到單元的擴(kuò)大剛度矩陣。單元?jiǎng)偠染仃囅禂?shù)取決于單元結(jié)點(diǎn)的局部編號(hào)順序，必須知道單元結(jié)點(diǎn)的局部編號(hào)與該結(jié)點(diǎn)在整體結(jié)構(gòu)中的總體編號(hào)之間的關(guān)系，才能得到單元?jiǎng)偠染仃囍械拿總€(gè)分塊在整體剛度矩陣中的位置。將單元?jiǎng)偠染仃囍械拿總€(gè)分塊按總體編碼順序重新排列后，可以得到單元的擴(kuò)大矩陣。假定單元結(jié)點(diǎn)的

19、局部編號(hào)與整體的對(duì)應(yīng)關(guān)系如下：?jiǎn)卧幪?hào)單元結(jié)點(diǎn)局部編號(hào)單元結(jié)點(diǎn)整體編號(hào)1i31j11m22i52j22m43i53j33m24i34j54m6單元（2）的單元擴(kuò)大矩陣的分塊矩陣形式如下，只列出非零的分塊：局部編號(hào)整體編號(hào)1234561234562）將全部單元的擴(kuò)大矩陣相加得到整體剛度矩陣。整體剛度矩陣如下所示：整體編號(hào)12345612+3+45+63.7約束條件的處理圖3.9所示的結(jié)構(gòu)的約束和載荷情況，如圖3.10所示。結(jié)點(diǎn)1、4上有水平方向的位移約束，結(jié)點(diǎn)4、6上有垂直方向的約束，結(jié)點(diǎn)3上作用有集中力（Px，Py）。圖3.10整體剛度矩陣K求出后，結(jié)構(gòu)上的結(jié)點(diǎn)力可以表示為：根據(jù)力的平衡，結(jié)點(diǎn)

20、上的結(jié)點(diǎn)力與結(jié)點(diǎn)載荷或約束反力平衡。用表示結(jié)點(diǎn)載荷和支桿反力，則可以得到結(jié)點(diǎn)的平衡方程：（3-29）這樣構(gòu)成的結(jié)點(diǎn)平衡方程組，在右邊向量P中存在未知量，因此在求解平衡方程之前，要根據(jù)結(jié)點(diǎn)的位移約束情況修改方程(3-29)。先考慮結(jié)點(diǎn)n有水平方向位移約束，與n結(jié)點(diǎn)水平方向?qū)?yīng)的平衡方程為：（3-30）根據(jù)支承情況，方程（3-30）應(yīng)該換成下面的方程：（3-31）對(duì)比公式（3-30）和（3-31），在式（3-29）中應(yīng)該做如下修正：在K矩陣中，第2n-1行的對(duì)角線元素改為1，該行中全部非對(duì)角線元素改為0；在P中，第2n-1個(gè)元素改為0。為了保持K矩陣的對(duì)稱性，將第2n-1列的全部非對(duì)角元素也改為0

21、。同理，如果結(jié)點(diǎn)n在垂直方向有位移約束，則（3-29）中的第2n個(gè)方程修改為， 在K矩陣中，第2n行的對(duì)角線元素改為1，該行中全部非對(duì)角線元素改為0；在P中，第2n個(gè)元素改為0。為了保持K矩陣的對(duì)稱性，將第2n列的全部非對(duì)角元素也改為0。（3-32）對(duì)圖3.9所示結(jié)構(gòu)的整體剛度在修改后可以得到以下的形式，（3-33）如果結(jié)點(diǎn)n處存在一個(gè)已知非零的水平方向位移，這時(shí)的約束條件為，（3-34）在K矩陣中，第2n-1行的對(duì)角線元素乘上一個(gè)大數(shù)A，向量P中的對(duì)應(yīng)換成，其余的系數(shù)保持不變。方程改為，（3-35）A的取值要足夠大，例如取1010。只有這樣，方程（3-35）才能與方程（3-34）等價(jià)。如果結(jié)

22、點(diǎn)n處存在一個(gè)已知非零的垂直方向位移，這時(shí)的約束條件為，。也可以采用同樣的方法修改整體剛度矩陣。3.8整體剛度矩陣的特點(diǎn)與存儲(chǔ)方法用有限元方法分析復(fù)雜工程問題時(shí)，結(jié)點(diǎn)的數(shù)目比較多，整體剛度矩陣的階數(shù)通常也是很高的。那么，是否在進(jìn)行計(jì)算時(shí)要保存整體剛度矩陣的全部元素？能否根據(jù)整體剛度矩陣的特點(diǎn)提高計(jì)算效率？整體剛度矩陣具有以下幾個(gè)顯著的特點(diǎn)：對(duì)稱性，稀疏性，非零系數(shù)帶形分布。1） 對(duì)稱性由單元?jiǎng)偠染仃嚨膶?duì)稱性和整體剛度矩陣的集成規(guī)則，可知整體剛度矩陣必為對(duì)稱矩陣。利用對(duì)稱性，只保存整體矩陣上三角部分的系數(shù)即可。2） 稀疏性單元?jiǎng)偠染仃嚨亩鄶?shù)元素為零，非零元素的個(gè)數(shù)只占較小的部分。如圖3.11所示

23、的結(jié)構(gòu)，結(jié)點(diǎn)2只和通過單元聯(lián)接的1、3、4、5結(jié)點(diǎn)相關(guān)，結(jié)點(diǎn)5只和通過單元聯(lián)接的2、3、4、6、8、9結(jié)點(diǎn)相關(guān)。由單元?jiǎng)偠染仃嚨奈锢硪饬x和整體剛度矩陣的形成方式可知，相關(guān)結(jié)點(diǎn)2、3、4、6、8、9及結(jié)點(diǎn)5本身產(chǎn)生位移時(shí)，才使結(jié)點(diǎn)5產(chǎn)生結(jié)點(diǎn)力，其余結(jié)點(diǎn)產(chǎn)生位移時(shí)不在該結(jié)點(diǎn)處引起結(jié)點(diǎn)力。在用分塊形式表示的整體矩陣中，與相關(guān)結(jié)點(diǎn)對(duì)應(yīng)的分塊矩陣具有非零的元素，其它位置上的分塊矩陣的元素為零，如圖3.12所示。圖3.11圖3.12整體剛度矩陣的分塊矩陣示意3）非零元素帶形分布在圖3.12中，明顯可以看出，整體剛度矩陣的非零元素分布在以對(duì)角線為中心的帶形區(qū)域內(nèi)，這種矩陣稱為帶形矩陣。在包括對(duì)角線元素的半個(gè)

24、帶形區(qū)域內(nèi)，每行具有的元素個(gè)數(shù)叫做半帶寬，用hbd表示。圖3.11所示結(jié)構(gòu)的相鄰結(jié)點(diǎn)編碼的最大差值為4，所以半帶寬為10。二維等帶寬存儲(chǔ)設(shè)整體剛度矩陣K為一個(gè)的矩陣，最大半帶寬為d。利用帶形矩陣的特點(diǎn)和對(duì)稱性，只需要保存以d為固定帶寬的上半帶的元素，稱為二維等帶寬存儲(chǔ)。進(jìn)行存儲(chǔ)時(shí)，把整體剛度矩陣K每行中的上半帶元素取出，保存在另一個(gè)矩陣K*的對(duì)應(yīng)行中，得到一個(gè)矩陣K*。把元素在K矩陣中的行、列編碼記為r、s，在矩陣K*中的行、列編碼記為r*、s*，對(duì)應(yīng)關(guān)系如下：r*=rs*=s-r+1 圖3.13(a)圖3.13(b)如圖3.13(a)所示的最大半帶寬為d的整體剛度矩陣K，采用二維等帶寬存儲(chǔ)后

25、得到如圖3.13(b)所示的矩陣K*。用新的方法存儲(chǔ)后，K矩陣中的對(duì)角線元素保存在新矩陣中的第1列中，K矩陣中的r行元素仍然保存在新矩陣的r行中，K矩陣中的s列元素則按照新的列編碼保存在新矩陣的不同列中。采用二維等帶寬存儲(chǔ)，需要保存的元素?cái)?shù)量與K矩陣中的總元素?cái)?shù)量之比為。所存儲(chǔ)的元素?cái)?shù)量取決于最大半帶寬d的值，d的值則由單元結(jié)點(diǎn)的編碼方式?jīng)Q定。雖然在采用二維等帶寬存儲(chǔ)時(shí)，仍然會(huì)保存一些零元素，但是采用這種方法時(shí)元素尋址很方便。圖3.14(a)圖3.14(b)對(duì)于同樣的有限元單元網(wǎng)格，按照?qǐng)D3.14(a)的結(jié)點(diǎn)編碼，最大的半帶寬為14；按照?qǐng)D3.14(b)的結(jié)點(diǎn)編碼，最大的半帶寬為18；按照?qǐng)D3

26、.11的結(jié)點(diǎn)編碼，最大的半帶寬為10。3.9線性方程組解法由于有限元分析需要使用較多的單元，線性方程組的階數(shù)很高，有限元求解的效率很大程度上取決于線性方程組的解法。利用矩陣的對(duì)稱、稀疏、帶狀分布等特點(diǎn)提高方程求解效率是關(guān)鍵。線性方程組的解法分為兩大類：直接解法，迭代解法。直接解法以高斯消去法為基礎(chǔ)，包括高斯消去法、等帶寬高斯消去法、三角分解法，以及適用于大型方程組求解的分塊算法和波前法等。迭代算法有高斯-賽德爾迭代、超松弛迭代和共軛梯度法等。在方程組的階數(shù)不是特別高時(shí)，通常采用直接解法。當(dāng)方程組的階數(shù)過高時(shí)，為避免舍入誤差和消元時(shí)有效數(shù)損失等對(duì)計(jì)算精度的影響，可以選擇迭代方法。這里給出了用Fo

27、rtran語(yǔ)言編寫的等帶寬高斯消去法的代碼，其中NROW為矩陣行的數(shù)目，NHBW為最大半帶寬。SUBROUTINE SOLVERB(NROW,NHBW,STIFF,DISPL)C Band elimination methodC .C : SOLVE EQUATIONS KS*H=(Q) :C :.:C KS is stored with the half band width methodC Input:CNROW - the quantity of rowsCNHBW -half band widthCSTIFFNROW,NHBWC - coefficient matrix stored with the half width methodCDISPLNROW - the right hand vectorC Output:CDISPL - results of variablesC Warning:Carray STIFF is changed.C*ADD:DPR* IMPLICIT DOUBLE PRECISION ( A-H,O-Z )C*END