桿系結(jié)構(gòu)有限元(4-5)

第1章桿系結(jié)構(gòu)有限元有限單元法1第1章概論1.1有限單元法及其發(fā)展概況1.2

有限元方法及解題步驟1.3

單元分析

1.4

整體分析1.5按單元定位向量形成總剛度方程1.6約束處理及求解2整體分析的任務(wù)：建立整體剛度方程[K]{D}={P}，主要目的是求解結(jié)點(diǎn)位移分量{D}。有限元位移法是在每一個(gè)結(jié)點(diǎn)上建立平衡方程，集合各結(jié)點(diǎn)的平衡方程得到一個(gè)平衡方程組[K]{D}={P}，出現(xiàn)在方程組內(nèi)的待定未知數(shù)便是求解的結(jié)點(diǎn)位移分量。31.4.1坐標(biāo)轉(zhuǎn)換矩陣在整體坐標(biāo)系中單元結(jié)點(diǎn)力向量和結(jié)點(diǎn)位移列向量可分別表示成坐標(biāo)變換關(guān)系為：4對(duì)于梁?jiǎn)卧鐖D所示，則有可簡(jiǎn)寫(xiě)為同理[T]稱(chēng)為坐標(biāo)變換矩陣5坐標(biāo)變換矩陣[T]坐標(biāo)變換矩陣的特點(diǎn)：從坐標(biāo)轉(zhuǎn)換矩陣[T]的表達(dá)式可以看出，[T]為正交矩陣，其逆矩陣等于其轉(zhuǎn)置矩陣，即有：

[T]-1=[T]T或者[T]T[T]=I

式中I為單位矩陣。于是有：6二、整體坐標(biāo)系下的單元?jiǎng)偠染仃囈胱鴺?biāo)變換關(guān)系，可得單元局部坐標(biāo)下的單元?jiǎng)偠确匠虨樽笥彝瑫r(shí)左乘[T]，即化簡(jiǎn)后可得整體坐標(biāo)下的剛度方程其中━━整體坐標(biāo)下的單元?jiǎng)偠染仃?對(duì)桿系結(jié)構(gòu)進(jìn)行單元分析，僅僅是有限元分析中的第一步。我們的目的是要對(duì)整個(gè)結(jié)構(gòu)進(jìn)行分析，研究結(jié)構(gòu)的整體性能。因此，在對(duì)結(jié)構(gòu)的各單元分析完成后，必須將單元分析的結(jié)果進(jìn)行整合，對(duì)結(jié)構(gòu)進(jìn)行整體分析。整體分析的過(guò)程實(shí)際上是如何將單元分析的結(jié)果進(jìn)行有效組合，建立整體剛度方程并求解結(jié)點(diǎn)位移的過(guò)程。根據(jù)對(duì)結(jié)點(diǎn)位移的編碼方式，可以采用“先處理法”和“后處理法”來(lái)建立整體剛度方程。1.4.2后處理法與先處理法8所謂后處理法，就是由單元?jiǎng)偠染仃囆纬烧w剛度矩陣，建立剛度方程后再引入支承條件，進(jìn)而求解結(jié)點(diǎn)位移的方法。運(yùn)用這種方法時(shí)，假設(shè)所有結(jié)點(diǎn)位移均為未知量，按照順序統(tǒng)一進(jìn)行編碼。后處理法9一、單元?jiǎng)偠确匠痰姆謮K表示

單元?jiǎng)偠确匠贪唇Y(jié)點(diǎn)編號(hào)分塊表示展開(kāi)寫(xiě)成:(1-4-1)(1-4-2)10二、總剛度方程的建立手算方法：根據(jù)桿端力和結(jié)點(diǎn)外載對(duì)每個(gè)結(jié)點(diǎn)建立平衡方程，如結(jié)點(diǎn)２的平衡條件可寫(xiě)為：外力單元①j端的桿端力單元②i端的桿端力11將式(1-4-2)的分塊形式代入F2①和F2②得（桿端力用結(jié)點(diǎn)位移表示）單元①j端的桿端力與1，2節(jié)點(diǎn)位移相關(guān)單元②i端的桿端力與2，3節(jié)點(diǎn)位移相關(guān)根據(jù)桿端位移與結(jié)點(diǎn)位移之間的諧調(diào)關(guān)系──代入幾何條件則12同理，分別取結(jié)點(diǎn)1,3,4由平衡條件得:13寫(xiě)成矩陣形式由此可見(jiàn)，總剛中的每一行（以結(jié)點(diǎn)分塊）的元素（子塊）是由與該結(jié)點(diǎn)相關(guān)單元的單剛元素集合（貢獻(xiàn)）而成的。結(jié)點(diǎn)123414程序方法：其實(shí)上述是一種手算方法，在程序設(shè)計(jì)中，并不象上面所述一樣，通過(guò)取每個(gè)結(jié)點(diǎn)列其平衡方程來(lái)得到總剛度方程，而是由單元結(jié)點(diǎn)編號(hào)確定每個(gè)單剛子塊在總剛中的位置──對(duì)號(hào)入座

15由此可見(jiàn)，每個(gè)單剛子塊元素在總剛中的位置，只需將單剛兩端的結(jié)點(diǎn)碼i,j（局部碼）換成結(jié)構(gòu)的結(jié)點(diǎn)編號(hào)碼（總碼），然后進(jìn)入總剛中的相應(yīng)位置即可，這一步驟又稱(chēng)為“換碼”，按此形成總剛的過(guò)程叫做“對(duì)號(hào)入座”形成總剛，這種方法又稱(chēng)直接剛度法。子塊元素的物理意義：Kij──j端發(fā)生單位位移時(shí)，i端上的抗力(剛度系數(shù))16例4-1：形成圖示剛架的總剛17解：?jiǎn)卧徒Y(jié)點(diǎn)編號(hào)如圖，由此知，總剛子塊由7×7塊組成。按對(duì)號(hào)入座方法得:18練習(xí)：按單剛節(jié)點(diǎn)分塊的形式，形成圖示結(jié)構(gòu)總剛。19先處理法利用先處理法對(duì)單元結(jié)點(diǎn)位移編碼時(shí)，僅對(duì)獨(dú)立的位移分量按自然數(shù)順序編號(hào)，若某些位移分量由于連接條件的限制彼此相等，則編為同一位移號(hào)，在支座處，由于剛性約束而使位移分量為零時(shí)，則對(duì)應(yīng)的編號(hào)為0。首先求出各單元的貢獻(xiàn)矩陣然后將它們疊加起來(lái)形成整體剛度矩陣202122232425第1章概論1.1有限單元法及其發(fā)展概況1.2

有限元方法及解題步驟1.3

單元分析

1.4

整體分析1.5按單元定位向量形成總剛度方程

1.6約束處理及求解261.5按單元定位向量形成總剛度方程

按單元定位向量形成總剛度方程前面介紹“對(duì)號(hào)入座”形成總剛的方法，是講子塊的對(duì)號(hào)入座，而在計(jì)算機(jī)程序中必須是將單剛的每個(gè)元素，用賦值語(yǔ)句送給總剛的相應(yīng)位置，這比子塊對(duì)號(hào)入座復(fù)雜，加上結(jié)構(gòu)各種不同的約束情況，使其更難處理。因此，在先處理法中，常引進(jìn)單元定位向量的概念。利用單元定位向量則可靈活地處理各種約束情況。單元定位向量就是確定單元各結(jié)點(diǎn)上位移分量編號(hào)的指示向量。形成單元定位向量的程序?qū)崿F(xiàn)方法很多，這里介紹一種利用結(jié)點(diǎn)定位向量形成單元定位向量的方法。27一、結(jié)點(diǎn)定位向量按結(jié)構(gòu)的結(jié)點(diǎn)編號(hào)順序，考慮各結(jié)點(diǎn)的約束信息和主從關(guān)系后，對(duì)每個(gè)結(jié)點(diǎn)所生成的位移分量編號(hào)叫結(jié)點(diǎn)定位向量。對(duì)于如圖所示結(jié)構(gòu)，其結(jié)點(diǎn)定位向量如圖中方括號(hào)所示(0表示有支承約束)。2829結(jié)點(diǎn)定位向量在程序中可用一個(gè)二維數(shù)組來(lái)儲(chǔ)存，如上圖，當(dāng)不儲(chǔ)存支座結(jié)點(diǎn)（三個(gè)位移分量均被約束）時(shí)為：（按結(jié)點(diǎn)編號(hào)順序存放）對(duì)于每個(gè)結(jié)點(diǎn)均由剛結(jié)點(diǎn)組成的剛架，當(dāng)考慮桿件的軸向變形時(shí)，每個(gè)結(jié)點(diǎn)３個(gè)位移分量，很容易由計(jì)算機(jī)自動(dòng)生成結(jié)點(diǎn)定位向量，結(jié)點(diǎn)i的第j個(gè)位移分量編號(hào)R為：K＝3*(i-1)+j30對(duì)于一般情況，在由計(jì)算機(jī)形成結(jié)點(diǎn)定位向量前必須告訴計(jì)算機(jī)哪些結(jié)點(diǎn)不是完全的剛結(jié)點(diǎn)──特殊結(jié)點(diǎn)，對(duì)每個(gè)特殊結(jié)點(diǎn)，其約束和主從關(guān)系都須告訴計(jì)算機(jī)，常將這種約束和主從關(guān)系信息叫做從結(jié)點(diǎn)信息（或特殊結(jié)點(diǎn)信息），將它們存放在從結(jié)點(diǎn)信息數(shù)組中：31約束信息可作如下規(guī)定：（1）獨(dú)立位移分量填0（2）支座約束填-1（3）位移與主結(jié)點(diǎn)相同填主結(jié)點(diǎn)號(hào)如上例共有７個(gè)從結(jié)點(diǎn)，其從結(jié)點(diǎn)信息為：從節(jié)點(diǎn)號(hào)X,Y,Z方向的約束信息32二、形成結(jié)點(diǎn)定位向量的子程序

INTEGERNTJ(90,4)COMMON/C17/JDW(150,3)OPEN(2,FILE='DATA')READ(2,*)NZ,N,NCJIF(NCJ.NE.0)READ(2,*)((NTJ(I,J),J=1,4),I=1,NCJ)CALLJDDWXL(NZ,N,NCJ,NTJ)WRITE(*,*)'WRITE(*,'(4I10)')(I,(JDW(I,J),J=1,3),I=1,NZ)END讀入從節(jié)點(diǎn)信息NZ：總節(jié)點(diǎn)數(shù)，N：可動(dòng)節(jié)點(diǎn)數(shù)，NCJ：特殊節(jié)點(diǎn)數(shù)33SUBROUTINEJDDWXL(NZ,N,NCJ,NTJ)COMMON/C17/JDW(150,3)DIMENSIONNTJ(90,4)DO10I=1,NZ’

NZ：總節(jié)點(diǎn)數(shù)DO10J=1,310JDW(I,J)=0’節(jié)點(diǎn)定位向量數(shù)組清零IF(NCJ.NE.0)THEN’假如特殊節(jié)點(diǎn)數(shù)不為零DO20I=1,NCJDO20J=1,320JDW(NTJ(I,1),J)=NTJ(I,J+1)ENDIFK=0DO40I=1,NDO40J=1,3從節(jié)點(diǎn)信息轉(zhuǎn)存到JDW數(shù)組34IF(JDW(I,J))50,60,70’小于零，等于零，大于零50JDW(I,J)=0‘有支承約束，不編號(hào)GOTO4060K=K+1‘該自由度可動(dòng)，編號(hào)加1JDW(I,J)=KGOTO40JDW(I,J)=JDW(JDW(I,J),J)‘從節(jié)點(diǎn)，信息>0，則找到對(duì)應(yīng)主節(jié)點(diǎn)的編號(hào)40CONTINUEEND第i節(jié)點(diǎn)是從節(jié)點(diǎn)，JDW(I,J)是從節(jié)點(diǎn)第J個(gè)位移分量對(duì)應(yīng)的主節(jié)點(diǎn)號(hào)，從主節(jié)點(diǎn)取出其第J個(gè)位移分量編號(hào)，賦給從節(jié)點(diǎn)。35前一頁(yè)的定位向量處理語(yǔ)句可改為(VB)：IFJDW(I,J)<0THENJDW(I,J)=0ELSEIFJDW(I,J)=0THENK=K+1JDW(I,J)=KELSEJDW(I,J)=JDW(JDW(I,J),J)ENDIF36三、單元定位向量

1.從單元信息數(shù)組中得到單元兩端的結(jié)點(diǎn)編號(hào)。（6）號(hào)單元：3，4節(jié)點(diǎn)2.由此編號(hào)便可從結(jié)點(diǎn)定位向量數(shù)組的對(duì)應(yīng)行取出即可構(gòu)成單元定位向量。（6）號(hào)單元的定位向量：

[MW]=[567809]可以只將結(jié)點(diǎn)定位向量一次生成并存入數(shù)組中，而單元定位向量則在需要時(shí)臨時(shí)構(gòu)成，以節(jié)約內(nèi)存。37形成單元定位向量子程序SUBROUTINEDWXL(E,MW)COMMON/C3/IJ(300,2)/C17/JDW(150,3)INTEGERE,MW(6)DO5I=1,65MW(I)=0I=IJ(E,1)DO10K=1,310MW(K)=JDW(I,K)J=IJ(E,2)DO20K=4,620MW(K)=JDW(J,K-3)ENDE：?jiǎn)卧?hào)單元E的i節(jié)點(diǎn)號(hào)單元E的j節(jié)點(diǎn)號(hào)38四、按單元定位向量形成總剛

由圖知，共有30個(gè)獨(dú)立位移分量，即總剛為30×30階，如取單元（6）為例，單剛元素入總剛的步驟為：1、形成單剛（6）2、形成單元定位向量[MW]（6）3、對(duì)單剛行、列循環(huán)，將其元素按定位向量指示的下標(biāo)對(duì)號(hào)入座進(jìn)總剛。39形成總剛子程序（半帶寬存儲(chǔ)）SUBROUTINEZG(M,U,L2,N)REAL*8A(450,66)COMMON/C5/A/C10/DK(6,6)/C17/JDW(150,3)INTEGERU,MW(6),EDO10I=1,MAXO(JDW(N,1),JDW(N,2),JDW(N,3)DO10J=1,ld10A(I,J)=0DO30E=1,M求最大自由度數(shù)N：可動(dòng)結(jié)點(diǎn)數(shù)40C調(diào)單剛CALLDG(E,12)C調(diào)單元定位向量25CALLDWXL(E,MW)C單剛元素入總剛DO30I=1,6DO30J=1,6IF(MW(I).LT.MW(J))GOTO30‘在總剛中列>行,即處于總剛的上三角部分,則丟棄IF(MW(I).EQ.0.OR.MW(J).EQ.0)GOTO30’被約束A(MW(I),MW(I)-MW(J)+1)=A(MW(I),MW(I)-MW(J)+1)+DK(I,J)30CONTINUEEND列號(hào)計(jì)算41半帶寬存儲(chǔ)（6）號(hào)單元42結(jié)構(gòu)總體剛度矩陣的特點(diǎn)1）對(duì)稱(chēng)性單剛對(duì)稱(chēng)，坐標(biāo)變換后單剛?cè)詫?duì)稱(chēng)2）奇異性后處理法集成的總剛具有奇異性3）稀疏帶狀性稀疏矩陣：大多數(shù)元素為零帶狀矩陣：只在對(duì)角線(xiàn)兩側(cè)一個(gè)狹窄的帶狀區(qū)域內(nèi)存在非零元素。43總剛的儲(chǔ)存方式在有限元法中，總剛總是一對(duì)稱(chēng)方陣，對(duì)于一個(gè)結(jié)點(diǎn)很多的結(jié)構(gòu)，其總剛的階數(shù)也很高，但往往是高度稀疏的，即0元素的比例很大，為了節(jié)省內(nèi)存，程序中一般都不以滿(mǎn)陣（方陣）儲(chǔ)存，常用的儲(chǔ)存方式主要有等帶寬半帶儲(chǔ)存（《實(shí)用有限單元法導(dǎo)論》稱(chēng)可省70％）和變帶寬一維儲(chǔ)存44最大半帶寬零元素零元素非零元素最大半帶寬主對(duì)角線(xiàn)元素等半帶存儲(chǔ)零元素存儲(chǔ)空間：nxn存儲(chǔ)空間：nxb45作業(yè)2：１、上述布置的自學(xué)任務(wù)２、將形成結(jié)點(diǎn)定位向量子程序加I/O(輸入/輸出)子程序要求：加輸入?yún)?shù)NA，用來(lái)控制以下幾種情況NA=0計(jì)梁柱軸向變形NA=1不計(jì)梁的軸向變形NA=2不計(jì)柱軸向變形NA=3不計(jì)梁柱軸向變形46練習(xí)：設(shè)用直接剛度法分析圖1所示平面框架，每個(gè)結(jié)點(diǎn)考慮水平、豎向、轉(zhuǎn)角3個(gè)位移分量(位移分量的方向如圖所示)，用先處理法形成的總剛度矩陣階數(shù)為(1)階，用后處理法形成總剛度矩陣階數(shù)為(2)階？先處理法圖2、圖3總剛階數(shù)？47若按圖示結(jié)點(diǎn)順序進(jìn)行自由度編號(hào)，每個(gè)節(jié)點(diǎn)的位移分量按水平、豎向、轉(zhuǎn)角位移順序排列，采用等半帶儲(chǔ)存總剛度矩陣，其半帶寬是(3)；設(shè)從下到上各樓層的層高分別為h1、h2、h3，從左到右