用有限元素法建立結構振動的數(shù)學模型

1、第三章 用有限元素法建立結構振動的數(shù)學模型31 引言【工程要求】：對于簡單的連續(xù)結構，如單件的桿、板、梁，可以建立結構振動的偏微分方程，但對于桿、板、梁組成的復雜結構，仍然采用建立偏微分方程的方法則十分困難。如果用假設模態(tài)法（李茲方法），對實際工程結構假設出品質(zhì)良好的整個結構的假設模態(tài)也十分困難。要對結構振動進行數(shù)值分析，必須建立振動的數(shù)學模型振動方程。工程結構振動分析中，要采用將結構離散為有限自由度系統(tǒng)的方法有限元素法，來建立結構的數(shù)學模型?！景l(fā)展簡況】有限元素法，是在上一世紀五十年代中期，經(jīng)過M.T.Turner及J.H.Argyris等人的開拓性工作以及后來許多研究者的大量工作，發(fā)展起來

2、的一種結構分析的有效方法，上一世紀六十年代初，由J.S.Archer及J.H.Argyris等人引入到結構動力學分析中來。有限元素法發(fā)展到今天，已經(jīng)非常成熟，而且與先進的計算機技術結合，已經(jīng)形成了一個以有限元分析方法為基礎的計算機輔助工程（CAE）的技術領域以及更進一步的虛擬產(chǎn)品設計（VPD）這樣的先進概念。世界上著名的CAE分析軟件商主要有MSC.software和Ansys等公司的產(chǎn)品?！居邢拊獎恿W分析的任務】在結構振動分析領域，有限元素法處理的問題主要是兩類：結構固有振動特性計算和結構振動響應計算（包括頻率響應分析與響應時間歷程分析）。兩類問題中，用有限元法建立振動數(shù)學模型是最基礎的工

3、作?！居邢拊胤ǎǚ治鼋Y構振動問題）的特點】：原則上，有限元素法由于其對復雜邊界的適應性，它可以處理任何復雜的結構。求解結果的精度可以根據(jù)需要不斷改善，建模過程規(guī)范統(tǒng)一，計算形式適合于計算機求解?！敬嬖诘膯栴}】：隨著精度要求的不斷提高，所要求的計算機容量和計算時間急劇增加，從而引出了大型特征值問題的快速求解方法、將大型結構振動問題轉化為若干小型結構振動問題集合的子結構求解方法，以及結構振動問題的并行求解方法等問題的研究?！竟こ探Y構振動分析方法】從結構振動分析的發(fā)展歷史看，經(jīng)典的方法有：1集中質(zhì)量法將質(zhì)量分別集中在若干節(jié)點處，形成集聚質(zhì)量陣。結構的剛度仍然連續(xù)分布，采用材料力學中求柔度的方法，求

4、出柔度系數(shù)，得到柔度矩陣，即用柔度法來形成剛度矩陣。集中質(zhì)量法存在問題：對大型復雜結構，用材料力學的方法，進行柔度矩陣的求解顯然是不現(xiàn)實的。2假設模態(tài)法以李茲法為基礎，選擇一組假設模態(tài)組成的模態(tài)矩陣，對結構進行離散，如第二章所述的方法。假設模態(tài)法存在的問題：1對幾何形狀復雜的結構，假設模態(tài)難以選擇。2對整個系統(tǒng)用假設模態(tài)法得到的運動方程是高度耦合的，求解困難。3對不同的結構，要根據(jù)實際情況選取不同的假設模態(tài)，求解過程不規(guī)范統(tǒng)一。引入有限元素法的思想既解決了上述方法的缺點，又保留了它們的優(yōu)點。【有限元法分析振動問題的基本原理】用有限元法分析結構動態(tài)問題的基本思想，與結構靜態(tài)分析的思想是一樣的。它

5、采用的方法仍然是：將結構分解為有限數(shù)目的單元，各元素間由節(jié)點相連，各單元內(nèi)結構的變形用位移形函數(shù)（相當于元素級的假設模態(tài)）來表示，以節(jié)點位移作為控制變量（元素的廣義坐標）。元素間的位移連續(xù)條件通過引入的形函數(shù)來滿足，動態(tài)平衡條件通過最后導出的有限元方程來體現(xiàn)。由于節(jié)點數(shù)目是有限的，最后得到的方程是一個多自由度、離散的、線性的矩陣微分方程。32 運動方程的建立仍然采用熟悉的拉格朗日方程法建立其數(shù)學模型（運動方程）。對任一單元內(nèi)部任一點的位移與節(jié)點位移的關系：（31）N稱為假設的已知位移形函數(shù)（可以看成是單元的假設模態(tài)，一般仍采用靜態(tài)變形函數(shù)）顯然： （32）單元的動能：（33）稱為單元質(zhì)量矩陣，

6、質(zhì)量陣是對稱矩陣。整個結構的動能為： （34）是全結構的節(jié)點位移列陣，為全結構質(zhì)量陣，代表對整個結構各單元的組集（Assemble）。單元的應變向量:， （35）幾何矩陣，彈性矩陣，應力向量單元的勢能為：（36）全結構的勢能： （37） （38）作用在單元上的分布力的虛功： （39）單元節(jié)點力（廣義力）（310）全結構的外力虛功：（311） （312）【阻尼的處理】采用粘性阻尼假定：阻尼力與運動速度成正比，方向與速度相反。單元中分布阻尼的耗散函數(shù)（瑞利耗散函數(shù)）： （313）耗散力（即瑞利耗散力）與耗散函數(shù)的關系為： （314）全結構的耗散函數(shù)： （315） （316）將全結構的動能、勢能、耗

7、散函數(shù)和廣義力代入非保守系統(tǒng)的拉格朗日方程， （317）得到： （318）【幾個相關問題】：1進行實際結構的振動分析時，在各個單元的矩陣組集之前，還要對單元矩陣進行由單元的局部坐標系向結構的總體坐標系轉換。記局部坐標系下節(jié)點位移向量向總體坐標系下節(jié)點位移向量的轉換陣為，則坐標轉換關系為： （319）其中為對單元矩陣組集時“對號入座”的定位矩陣。（320）（321）（322） （323）經(jīng)過上述變換后的單元矩陣可以直接疊加得到結構總體矩陣。2為了求解結構的固有振動特性，需要求解無阻尼情況下結構的自由振動方程：(324)將固有振動的簡諧運動形式（325）代入得到的結構的特征方程：或 （326）數(shù)學

8、上構成所謂的廣義特征值問題。3振動分析中采用的質(zhì)量陣問題在結構振動分析中，常采用的質(zhì)量陣形成方法有：集中質(zhì)量模型和一致質(zhì)量模型。在采用集中質(zhì)量模型時，一般是按照杠桿原理將單元質(zhì)量向單元各個節(jié)點上進行分配，在局部轉動效應顯著時，還要考慮單元的轉動慣量。集中質(zhì)量模型得到的質(zhì)量陣為對角矩陣。將單元內(nèi)慣性分布視為與靜力形函數(shù)同樣規(guī)律的分布，導出的質(zhì)量矩陣稱為一致質(zhì)量陣。即上面（）推導出的質(zhì)量陣。這樣的質(zhì)量陣為滿陣。注意：集中質(zhì)量陣與一致質(zhì)量陣都不是振動結構在實際上精確的質(zhì)量分布模型。理論上，結構的動位移是與頻率相關的。動位移在不同振動頻率和振型下是不同的，即不同的頻率對應有不同的慣性。所以嚴格地講，質(zhì)

9、量陣也是與頻率相關的。下面以軸向振動的桿元為例，說明這個問題。作為連續(xù)體的二力桿元的振動偏微分方程（波動方程）為：xL12 （327）記稱為波速，分別為桿的彈性模量和密度。方程（327）的解為：（328）故此時形函數(shù)為: （329）代入單元質(zhì)量陣公式得：（330）顯然，這個質(zhì)量陣為滿陣，且各元素為與頻率相關的量。理論上這樣的質(zhì)量陣真實反映實際的慣性分布，計算得到的固有特性更精確，但卻使計算大大復雜化，而且只有對簡單的情況（如二力桿元）下，才能從偏微分振動方程得出（動力）形函數(shù)，因此，工程實際中很少采用這種方法。而是采用與推導剛度陣時一致的（靜力）形函數(shù)。至于采用集中質(zhì)量陣還是一致質(zhì)量陣，得到的

10、結果更好，沒有固定的規(guī)律和結論，要視具體情況而定。從一般經(jīng)驗上講，在單元劃分較細時，用集中質(zhì)量陣較好，反之宜采用一致質(zhì)量陣。根據(jù)特征值隔離定理知道，采用一致質(zhì)量陣分析得到的是固有頻率精確解的上界，隨著分元的細化，計算結果單調(diào)地向精確解逼近，而集中質(zhì)量陣給出的結果就不具備這種特性，可能偏低也可能偏高。從工程分析經(jīng)驗看，由于建立有限元模型的離散過程，已經(jīng)使結構比實際結構的剛度增大，因此采用集中質(zhì)量陣可能有時反而會得到誤差較小的結果，但這需要經(jīng)驗和技巧。但采用集中質(zhì)量陣得到的振型一般誤差較大，對振型要求較高時，還是宜采用一致質(zhì)量陣。集中質(zhì)量陣是對角陣，在計算時可以節(jié)省計算時間。我們希望能獲得一種優(yōu)于

11、集中質(zhì)量陣的對角化質(zhì)量陣。例如，對于梁的彎曲振動，可以按下式來計算對角化的非一致質(zhì)量陣的單元： （331）而非對角元全部置零，是梁的四個形函數(shù)。由于這種質(zhì)量矩陣在一定程度上反映了單元形函數(shù)的特性，因此可以給出精度較好的計算結果。33 典型結構單元的有限元建模結構有限元模型的自由度數(shù)節(jié)點數(shù)×節(jié)點位移數(shù)一、 縱向振動桿元桿元的任意截面處位移由兩節(jié)點位移插值得到。 （332）應滿足邊界條件： （333）桿元受節(jié)點力時的靜態(tài)方程為： （334）其解為： （335）代入邊界條件得到形函數(shù)： （336）從而桿元的質(zhì)量陣和剛度陣為： （337） （338）二、 彎曲振動梁元梁元的任意截面處位移由兩

12、節(jié)點位移插值得到。（339）為梁元的節(jié)點位移。均勻梁元受常值節(jié)點力作用時的撓曲線偏微分方程為： （340）從而 （341）應滿足邊界條件： （342）代入邊界條件得到： （343）代入位移表達式，整理后，得到形函數(shù)為： （344）形函數(shù)陣： （345）梁元的質(zhì)量陣和剛度陣為： （346） （447）上面給出的是平面梁元的特性矩陣，對于空間梁元特性矩陣形成過程完全相同?？梢詤⒖既魏我槐居邢拊胤ǚ矫娴膶Ｖ蚪滩摹＿@里不再贅述?！締卧淖鴺俗儞Q陣】：單元局部坐標與總體坐標不一致時，需要進行坐標變換，記局部坐標與總體坐標之間夾角為，則局部坐標下節(jié)點位移列陣與總體坐標下位移列陣間的變換關系為： (34

13、8)則坐標變換陣為： (349)三面內(nèi)振動的平板元以三節(jié)點三角元為例，取三角元三個頂點的位移為單元廣義坐標， （350）單元內(nèi)任一點的位移為：（351）三個形函數(shù)應滿足邊界條件： （352）設滿足此條件的形函數(shù)為： （353）代入邊界條件可得到： （354） （355）為三角形單元的面積。由上式按下圖輪換下標求得：（356a）幾何矩陣：（356） （357）對平面應力問題： (358a) （358） （359）質(zhì)量陣和剛度矩陣為： （360）記板厚為，且由于都是常數(shù)矩陣，則 （361）的分塊形式為 （362）節(jié)點剛度陣（363）單元質(zhì)量陣（364）四彎曲振動的平面板元仍采用三節(jié)點三角元，引入局

14、部坐標系，原點在單元的1節(jié)點，軸沿單元的一條邊，軸垂直于單元平面，如圖所示。單元的三個節(jié)點共有九個節(jié)點位移條件，故設單元內(nèi)任意一點的離面垂直位移表達式為： （365）或 （366）其中 （367）（368）單元的節(jié)點位移邊界條件為： （369） （365）代入（369）式得到： （370） （371）顯然，各節(jié)點的離面位移垂直位移（撓度）、繞軸的轉角、繞軸的轉角，就是單元特性分析時所用的廣義坐標。節(jié)點的位移列陣為： （372）從而有：， （373）因此，彎曲板元的形函數(shù)為： （374）按薄板彎曲問題的基本假定（直法線假定），板內(nèi)各點線位移為： （375）應變列陣： （376） （377） （3

15、78）其中 （379）為微分算子列陣。記（380）則 （381）從而得到彎曲板元的幾何矩陣：（382）單元剛度陣 （383）（384）為彈性矩陣。單元質(zhì)量矩陣為： （385）（具體的計算結果略）積分后單元的剛度矩陣和質(zhì)量矩陣表達式可參見有限元素法的專著。其它形式的單元的特性矩陣，均可按上述過程推導出。當然，在實際結構振動有限元分析中，如果自己編程，不需要自行推導，可以查閱有關參考書籍。應該承認，在各種關于結構有限元靜力分析書籍中，對有限元的介紹非常詳盡，也給出了各種高精度的單元形函數(shù)和剛度矩陣，但一般有限元素法的書籍中沒有給出單元質(zhì)量矩陣。如果采用商用有限元分析軟件，則在選擇單元類型時，要根據(jù)

16、具體結構、對分析結果的精度要求以及網(wǎng)格劃分的情況，來選擇最合適的單元。這些問題在有限元著作中都有詳細論述，這里也不再贅述。34結構振動分析有限元建模中的幾個問題【阻尼矩陣的工程處理】在對結構進行振動響應分析時，必須計及阻尼的影響，然而由于實際結構的阻尼特性比較復雜，即使我們可以借助瑞利粘性阻尼假設，得到帶粘性阻尼結構的振動方程（386）但是，結構的阻尼特性分布參數(shù)無法從理論上確定，因而就無法用理論方法來計算矩陣中的元素。在工程振動分析中，常常是采用瑞利比例阻尼假設，即假定阻尼陣是質(zhì)量陣和剛度陣的線性組合： （387）比例常數(shù)通過對結構或類似結構，進行模態(tài)試驗，測得其第階和第階模態(tài)頻率和阻尼比來

17、確定。理論上可以從下面方程解出：（388）一般是測定多階模態(tài)參數(shù)，通過最小二乘法求解出。【剛度矩陣奇異性的處理】由于采用有限單元法，通過對各單元剛度矩陣的裝配，得到的全結構總體剛度陣一般是無約束結構的剛度矩陣，它是奇異的，即不能求逆。因此需要按照結構實際的約束邊界條件進行處理，在對剛度矩陣進行處理時，有四種方法：1 第一種方法：當給定的邊界約束條件足以使結構成為靜定或超靜定結構時，只需在進行全結構的總剛度陣裝配時，抽掉剛度矩陣中對應于結構節(jié)點位移為零的行和列，即得到一個非奇異的降階的剛度矩陣，對質(zhì)量陣進行相應的處理。這樣處理的缺點是不能進一步求出約束反力。2 第二種方法：用無限大或計算機允許的

18、最大數(shù)來代替剛度陣中對應于零位移節(jié)點的主對角元素，從而使剛度陣成為非奇異的可逆矩陣，并使對應的節(jié)點位移具有一個非常小、十分接近于零的解。這樣處理的缺點是未能利用邊界節(jié)點約束降低結構矩陣階數(shù)而不能減小計算工作量。3 第三種方法：對剛度矩陣和質(zhì)量矩陣按節(jié)點位移是否為零進行分塊，如對結構的自由振動方程 （389）為約束反力列陣，下表，分別代表自由節(jié)點和約束節(jié)點。將矩陣方程分塊展開得到： （390）解第一式的特征值問題得到及，再代入第二式得到約束反力。顯然由這樣的方法可以得到結構的全部信息，但是增加了編程計算的難度。4 當結構屬于自由自由結構，或結構具有剛體運動模態(tài)時，我們不可能根據(jù)結構的邊界約束條件

19、來消除剛度矩陣的奇異性，從而給特征值問題的求解帶來困難。為了消除剛度矩陣的奇異性，可以利用特征值問題的移軸特性，將剛度矩陣轉換成正定的矩陣。如下所示： （391）為一個任意的、適當?shù)恼?shù)，上方程中是一個正定矩陣，可以按通常的方法求解，得到的頻率與原來結構的固有頻率間具有關系： （392）而振型不變，就是由上方程解得的。值得注意的是，由于結構具有個剛體模態(tài)，對應的固有頻率為零，所以，上方程的解也包括個剛體模態(tài)，對應著個重頻（這里用了二階小量符號，因為一般數(shù)值求解的結果不一定得到理論上的零），在求解時可以用“掃模法”將這個剛體模態(tài)掃去后再用一般的求解特征值方法（程序）求解?！咀杂啥鹊撵o力縮聚法】為了提高計算精度，在結構的有限元建模時，對結構可以采用網(wǎng)格細分的方法得到較精細的有限元網(wǎng)格；也可以采用較粗的網(wǎng)格劃分，而采用帶有內(nèi)節(jié)點的高階元素。在采用這種高階元素時，內(nèi)節(jié)點參數(shù)與單元間的位移連續(xù)條件是無關的，與單元本身以外的節(jié)點參數(shù)也沒有直接聯(lián)系。因此可以用某種方法在組成全結構矩陣前，將其消去以減少結構的總自由度數(shù)。還有一種情況，如對梁彎曲振動問題進行分析時，往往只保留節(jié)點橫向位移參數(shù)，而略去節(jié)點的轉角參數(shù)，也需要在組成全結構的矩陣后，用節(jié)點橫向位移參數(shù)來表示節(jié)點轉角參數(shù)，從而達到特征值問題降階的目的。此外，許