用于振動(dòng)分析的有限元方法課件

用于振動(dòng)分析的有限元方法指導(dǎo)老師：陳益報(bào)告人：成志斌韓宗彪何瑜寧鵬用于振動(dòng)分析的有限元方法指導(dǎo)老師：陳益內(nèi)容有限元介紹單個(gè)元素的運(yùn)動(dòng)方程整個(gè)系統(tǒng)的運(yùn)動(dòng)方程整個(gè)系統(tǒng)的邊界條件的加載及質(zhì)量矩陣

MATLAB實(shí)例及總結(jié)單個(gè)元素的質(zhì)量矩陣、剛度矩陣、力矢量及其轉(zhuǎn)化內(nèi)容有限元介紹單個(gè)元素的運(yùn)動(dòng)方程整個(gè)系統(tǒng)的運(yùn)動(dòng)方程整個(gè)系統(tǒng)的有限元法簡(jiǎn)介有限元法是一種可用于精確地(但近似)解決許多復(fù)雜的振動(dòng)問題的數(shù)值方法。對(duì)于基本的一維元素進(jìn)行有限元分析，能得到質(zhì)量矩陣與剛度矩陣和所需的力矢量，對(duì)于二維三維，元素矩陣會(huì)轉(zhuǎn)換成相關(guān)的更高維的空間。使用一致的和集中質(zhì)量矩陣的有限元方程并結(jié)合邊界條件能為復(fù)雜系統(tǒng)提供解釋。最后,使用MATLAB程序得到在軸向載荷下的指定節(jié)點(diǎn)位移,固有振動(dòng)頻率和特征值分析。有限元法簡(jiǎn)介有限元法是一種可用于精確地(但近似本章目的*認(rèn)識(shí)用于解決不同類型振動(dòng)問題的剛度和質(zhì)量矩陣。*將矩陣元素從局部坐標(biāo)系變換到全球坐標(biāo)系。*裝配單元矩陣和應(yīng)用邊界條件。*對(duì)桿、梁元素進(jìn)行靜態(tài)分析。*對(duì)桿、梁元素進(jìn)行動(dòng)態(tài)分析來得到固有頻率和振型。*在有限元振動(dòng)分析使用一致的集中質(zhì)量矩陣。*使用MATLAB解決振動(dòng)問題。本章目的*認(rèn)識(shí)用于解決不同類型振動(dòng)問題的剛度和質(zhì)量矩陣。有限元思想1，實(shí)際結(jié)構(gòu)被一些元素所取代,這些元素都是被假定為一個(gè)連續(xù)的結(jié)構(gòu)部件即有限元，這些元素在特定點(diǎn)即節(jié)點(diǎn)上互相關(guān)聯(lián)。2，如果解決方案的各方面都選擇得當(dāng),那么它可以收斂到精確的解決方案，因?yàn)榻M成總體結(jié)構(gòu)的元素很小，在節(jié)點(diǎn)上的力的平衡和元素之間的位移都令人感到滿意,這樣整個(gè)結(jié)構(gòu)(組合的元素)表現(xiàn)為單一實(shí)體。3，因?yàn)榈玫綔?zhǔn)確解很難，所以得到一個(gè)方便且逼近的近似解很有價(jià)值。有限元思想1，實(shí)際結(jié)構(gòu)被一些元素所取代,這些元素都是被假定為元素的運(yùn)動(dòng)方程龍門刨銑床有限元模型三角板元素梁元素元素的運(yùn)動(dòng)方程龍門刨銑床有限元模型三角板元素梁元素元素的運(yùn)動(dòng)方程位移函數(shù)形狀函數(shù)各點(diǎn)對(duì)應(yīng)位移未知節(jié)點(diǎn)位移數(shù)n動(dòng)能應(yīng)變能質(zhì)量矩陣剛度矩陣位移函數(shù)形狀函數(shù)各點(diǎn)對(duì)應(yīng)位移未知節(jié)點(diǎn)位移數(shù)n動(dòng)能應(yīng)變能元素的運(yùn)動(dòng)方程位移函數(shù)質(zhì)量矩陣剛度矩陣位移函數(shù)主要內(nèi)容：一，單元的質(zhì)量、剛度矩陣、等效節(jié)點(diǎn)力二，單元矩陣的坐標(biāo)變換三，整個(gè)系統(tǒng)的運(yùn)動(dòng)方程主要內(nèi)容：一單元的質(zhì)量、剛度矩陣，等效節(jié)點(diǎn)力矢量一桿單元一個(gè)桿單元是從桿上劃分出的一個(gè)小段，如下圖所示。由于單元很小，ρ、A均視為常量?，F(xiàn)在就以這最簡(jiǎn)單的桿單元，推導(dǎo)出它的質(zhì)量、剛度矩陣，等效節(jié)點(diǎn)力。圖12.1圖12.1一單元的質(zhì)量、剛度矩陣，等效節(jié)點(diǎn)力矢量一桿單元一個(gè)桿單元(1)求桿單元上任意點(diǎn)的位移u(x,t)本來，桿單元上任意點(diǎn)的位移u(x,t)與節(jié)點(diǎn)的位移u1(t)、u2(t)之間的關(guān)系是未知的，但是，只要單元?jiǎng)澐值淖銐蛐?，那么其間的關(guān)系就無關(guān)大局。所以可以假定它們之間有簡(jiǎn)單的線性關(guān)系，即根據(jù)節(jié)點(diǎn)位移對(duì)單元內(nèi)任意點(diǎn)位移進(jìn)行插值：（1）式中，Φ1、Φ2稱為線性系數(shù)，與單元里點(diǎn)的位置有關(guān)，是x的函數(shù)。此函數(shù)與單元的形狀有關(guān)，又叫形狀函數(shù)。(1)求桿單元上任意點(diǎn)的位移u(x,t)本來，桿單元上任意形狀函數(shù)和插值函數(shù)一樣是任意的，但必須邊界條件：只有滿足此條件單元才能協(xié)調(diào)一致運(yùn)動(dòng)，而不致破壞系統(tǒng)的完整性，因此這兩個(gè)條件實(shí)際上就是變形協(xié)調(diào)條件。將式（1）帶入（2）中，就可以得到形狀函數(shù)Φ1(x)、Φ2(x)所滿足的邊界條件：（2）（3）形狀函數(shù)和插值函數(shù)一樣是任意的，但必須邊界條件：只有滿足此條以上邊界條件確定了由于這兩個(gè)函數(shù)的任意性我們可以用簡(jiǎn)單的線性函數(shù)來近似，因此有：（4）代回（1）式中有：（5）為此，我們已經(jīng)找到了用節(jié)點(diǎn)位移表示單元內(nèi)任意一點(diǎn)位移的表達(dá)式。以上邊界條件確定了由于這兩個(gè)函數(shù)的任意性我們可以用簡(jiǎn)單的線性（2）計(jì)算此單元的動(dòng)能和勢(shì)能桿單元的動(dòng)能可表示成：（6）（2）計(jì)算此單元的動(dòng)能和勢(shì)能桿單元的動(dòng)能可表示成：（6）上式中，ρ是材料的密度，A是桿單元的橫截面積。用矩陣形式表示（6）式為：（7）其中，（8）所以，質(zhì)量矩陣可以認(rèn)為是：（9）上式中，ρ是材料的密度，A是桿單元的橫截面積。用矩陣形式表示桿單元的勢(shì)能可以寫成：（10）式中，E是彈性模量，（10）表示成矩陣形式為：桿單元的勢(shì)能可以寫成：（10）式中，E是彈性模量，（10）表這里，，所以剛度矩陣[k]可以表示成：（11）（12）(3)計(jì)算等效節(jié)點(diǎn)力設(shè)單元上x處作用有分布力f(x,t)，現(xiàn)在要把它等效成節(jié)點(diǎn)力遵循等效原則，即原載荷和等效之后的節(jié)點(diǎn)載荷在虛位移上所做的虛功相等。這里，，所以剛度矩陣[k]可以表示成：（11）（12）(3)其實(shí)，就是對(duì)應(yīng)于廣義坐標(biāo)的廣義力，為此，計(jì)算所做的虛功：把上式寫成矩陣形式：（13）其實(shí)，就是對(duì)應(yīng)于廣義坐標(biāo)的廣義力，為此，計(jì)算所做的虛功：把上所以等效節(jié)點(diǎn)力可以寫成：（14）所以等效節(jié)點(diǎn)力可以寫成：（14）二梁?jiǎn)卧缦聢D所示，一個(gè)梁?jiǎn)卧彩怯袃蓚€(gè)節(jié)點(diǎn)，但是有四個(gè)自由度，每個(gè)節(jié)點(diǎn)處，有兩種位移形式，一個(gè)是線位移，即撓度，一種是角位移。圖中，是力，是力矩。是分布載荷是對(duì)應(yīng)的線位移，是對(duì)應(yīng)的轉(zhuǎn)角。是梁?jiǎn)卧先我馕灰苮處的撓度。圖12.2二梁?jiǎn)卧缦聢D所示，一個(gè)梁?jiǎn)卧彩怯袃蓚€(gè)節(jié)點(diǎn)，但是有四個(gè)自在靜載彎曲條件下，梁?jiǎn)卧先我恻c(diǎn)出的撓度是x的三次方程，可寫成：此方程必須滿足下面的邊界條件：由此可以求解處a(t)、b(t)、c(t)、d(t)，進(jìn)而撓度方程為：（16）（17）（18）在靜載彎曲條件下，梁?jiǎn)卧先我恻c(diǎn)出的撓度是x的三次方程，可寫上式可以寫成形狀函數(shù)的表示：其中，形函數(shù)分別為：梁?jiǎn)卧膭?dòng)能、勢(shì)能、虛功表達(dá)式分別為：（19）上式可以寫成形狀函數(shù)的表示：其中，形函數(shù)分別為：梁?jiǎn)卧膭?dòng)能式中I是橫截面的慣性矩上式中：（20）（21）（22）式中I是橫截面的慣性矩上式中：（20）（21）（22）通過上式，可以得到梁?jiǎn)卧馁|(zhì)量、剛度矩陣，等效節(jié)點(diǎn)力：通過上式，可以得到梁?jiǎn)卧馁|(zhì)量、剛度矩陣，等效節(jié)點(diǎn)力：二單元矩陣的坐標(biāo)變換局部坐標(biāo)系：以各個(gè)單元本身的軸線為基準(zhǔn)所設(shè)立的坐標(biāo)系。便于計(jì)算節(jié)點(diǎn)位移。缺點(diǎn)：如果整個(gè)系統(tǒng)里各個(gè)單元取向各異，各個(gè)節(jié)點(diǎn)位移方向不一致，如下圖。如何使匯交于一個(gè)節(jié)點(diǎn)的各個(gè)桿件的節(jié)點(diǎn)位移真正相等?解決方法：進(jìn)行坐標(biāo)變換如右圖的系統(tǒng)，有四個(gè)桿件，u1(t)、u2(t)為局部坐標(biāo)系的節(jié)點(diǎn)位移，Ui

為全局坐標(biāo)系下的位移圖12.3二單元矩陣的坐標(biāo)變換局部坐標(biāo)系：以各個(gè)單元本身的軸線為基準(zhǔn)如下圖，節(jié)點(diǎn)位移在局部、全局坐標(biāo)系中的關(guān)系：坐標(biāo)變換矩陣（23）如下圖，節(jié)點(diǎn)位移在局部、全局坐標(biāo)系中的關(guān)系：坐標(biāo)變換矩陣（2其中，因?yàn)閱卧膭?dòng)能、勢(shì)能與坐標(biāo)系無關(guān):(24)其中，因?yàn)閱卧膭?dòng)能、勢(shì)能與坐標(biāo)系無關(guān):(24)得到在全局坐標(biāo)下的單元質(zhì)量、剛度矩陣為：類似地，根據(jù)單元在兩個(gè)坐標(biāo)系下的力所做的虛功相等：得到在全局坐標(biāo)系下的等效節(jié)點(diǎn)力：得到在全局坐標(biāo)下的單元質(zhì)量、剛度矩陣為：類似地，根據(jù)單元在兩三全系統(tǒng)運(yùn)動(dòng)方程經(jīng)過坐標(biāo)變換，各個(gè)單元的節(jié)點(diǎn)位移方向被統(tǒng)一起來，但是不同的節(jié)點(diǎn)有不同的節(jié)點(diǎn)位移，為了便于綜合出全系統(tǒng)的運(yùn)動(dòng)方程，首先要建立全系統(tǒng)的節(jié)點(diǎn)位移向量。每個(gè)單元的節(jié)點(diǎn)位移向量與全系統(tǒng)的節(jié)點(diǎn)位移向量之間的關(guān)系：長(zhǎng)方形矩陣由1和0組成單元節(jié)點(diǎn)位移向量單元節(jié)點(diǎn)位移向量三全系統(tǒng)運(yùn)動(dòng)方程經(jīng)過坐標(biāo)變換，各個(gè)單元的節(jié)

例如，圖的12.5中的1單元，方程，變?yōu)椋喊衙總€(gè)單元的動(dòng)能相加，就得到了整個(gè)系統(tǒng)的動(dòng)能：（把整個(gè)系統(tǒng)的動(dòng)能表示成關(guān)于節(jié)點(diǎn)速度矢量的形式）例如，圖的12.5中的1單元，方程，變?yōu)椋喊衙總€(gè)單元的動(dòng)能這樣就得到了整個(gè)系統(tǒng)的質(zhì)量矩陣:類似地，考慮整個(gè)系統(tǒng)的勢(shì)能，便可以得到整個(gè)系統(tǒng)的剛度矩陣：整個(gè)系統(tǒng)的廣義力向量：最后得到整個(gè)系統(tǒng)的運(yùn)動(dòng)方程：這樣就得到了整個(gè)系統(tǒng)的質(zhì)量矩陣:類似地，考慮整個(gè)系統(tǒng)的勢(shì)能，12.6——添加邊界條件前文中，節(jié)點(diǎn)沒有固定，結(jié)構(gòu)在節(jié)點(diǎn)力的作用下會(huì)發(fā)生剛體位移。也就是說，矩陣[K]是奇異矩陣。通常情況下，我們希望結(jié)構(gòu)的位移為零。

因此，我們需要添加邊界條件對(duì)矩陣[M]、[K]和向量F進(jìn)行約束。N：結(jié)構(gòu)中自由節(jié)點(diǎn)位移的數(shù)目12.6——添加邊界條件前文中，節(jié)點(diǎn)沒有固定，結(jié)構(gòu)例1：桿件分析如圖：均質(zhì)；長(zhǎng)0.5m；斷面截面積5e-4m^2；楊氏模量200GPa；密度7850Kg/m^3；左端固定。a.節(jié)點(diǎn)2處施加1000N靜態(tài)軸向外力u2，求應(yīng)力b.求系統(tǒng)固有頻率例1：桿件分析如圖：均質(zhì)；長(zhǎng)0.5m；斷面截解：a：平衡方程：A=5e-4，E=2e11，l=0.5，f2=1000，代入方程得：u1：位移，f1：節(jié)點(diǎn)1處應(yīng)力，添加邊界條件：u1=0，解得： u2=5e-10m解：由應(yīng)力與應(yīng)變的關(guān)系：

表示長(zhǎng)度的變化，表示應(yīng)變；b：

由剛度矩陣和質(zhì)量矩陣，得特征值方程：

由應(yīng)力與應(yīng)變的關(guān)

式中w位固有頻率，U1、U2分別是節(jié)點(diǎn)1、2的振幅，添加邊界條件：U1=0；解得：式中w位固有頻率，U1、U2分別是節(jié)點(diǎn)1

例2：梁的自然頻率解：梁被理想化為單一單元，局部和整體的節(jié)點(diǎn)位移相同，如圖所示：

梁的剛度矩陣：例2：梁的自然頻率解：梁被理想化為單一單元，局部

質(zhì)量矩陣：

節(jié)點(diǎn)位移向量：

與端點(diǎn)相關(guān)的邊界條件：W1=0，W3=0；解得：質(zhì)量矩陣：節(jié)點(diǎn)位移向量：

求解特征值：

乘以l/2EI得：

令系數(shù)矩陣的行列式等于0得：

方程的根即梁的自然頻率：求解特征值：乘以l/2EI得

結(jié)果可以和精確解比較：結(jié)果可以和精確解比較：12.7——一致、集中質(zhì)量矩陣

12.3節(jié)中推出的質(zhì)量矩陣是一致質(zhì)量矩陣，因?yàn)橛糜谕茖?dǎo)剛度矩陣的位移模型也用于推導(dǎo)質(zhì)量矩陣。

一些動(dòng)態(tài)問題可以用形式簡(jiǎn)單的質(zhì)量矩陣求解。最簡(jiǎn)單的質(zhì)量矩陣——集中質(zhì)量矩陣，可以通過將質(zhì)點(diǎn)指定到節(jié)點(diǎn)上。

集中質(zhì)量針對(duì)平移和旋轉(zhuǎn)的元素,假設(shè)在平均位置兩側(cè)的特定位移表現(xiàn)得像個(gè)剛體而剩余的元素不參與運(yùn)動(dòng)。

因此這種假設(shè)不包括元素位移之間存在的動(dòng)態(tài)耦合,因此產(chǎn)生的元素質(zhì)量矩陣是純粹的對(duì)角矩陣。12.7——一致、集中質(zhì)量矩陣12.3節(jié)中1、桿的集中質(zhì)量矩陣：2、梁的集中質(zhì)量矩陣：

旋轉(zhuǎn)自由度的慣性影響被假定為0；若考慮慣性影響，有轉(zhuǎn)動(dòng)慣量：

集中質(zhì)量矩陣變?yōu)椋?/p>

1、桿的集中質(zhì)量矩陣：

對(duì)于一般的動(dòng)態(tài)問題，兩者誰能得到更精確的解？

兩個(gè)質(zhì)量矩陣很相似，他們不考慮各種位移自由度的元素之間的動(dòng)態(tài)耦合。

他們的形狀函數(shù)也近似，都是用靜態(tài)的位移模型推導(dǎo)而來。

然而,由于集中質(zhì)量矩陣對(duì)角，在計(jì)算時(shí)他使用更少的存儲(chǔ)空間。

下面的例子說明了在一個(gè)簡(jiǎn)單的振動(dòng)問題中，集中和一致質(zhì)量矩陣的應(yīng)用。集中質(zhì)量矩陣與一致質(zhì)量矩陣：對(duì)于一般的動(dòng)態(tài)問題，兩者誰能得到更精確的解例：桿的一致和集中質(zhì)量矩陣

用一致和集中質(zhì)量矩陣求如圖所示兩端固定桿的固有頻率，用兩個(gè)桿單元建模。解：?jiǎn)卧膭偠群唾|(zhì)量矩陣分別是：

質(zhì)量矩陣的下標(biāo)c和l分別表示一致和

集中質(zhì)量矩陣。例：桿的一致和集中質(zhì)量矩陣用一致由于該桿由兩個(gè)單元建模，組合的剛度和質(zhì)量矩陣如下：

方框中的部分分別與單元1和2相關(guān)。

由于該桿由兩個(gè)單元建模，組合的剛度和質(zhì)量矩陣添加邊界條件U1=U3=0后，特征值問題為：特征值w由以下方程課解：代入已知條件得：

用一致質(zhì)量矩陣

用集中質(zhì)量矩陣

添加邊界條件U1=U3=0后，特征值問題為：

事實(shí)上，方程的精確解析解是：解得：事實(shí)上，方程的精確解析解是：12.8MATLAB應(yīng)用舉例例12.5階梯軸的有限元分析圖12.11中的階梯軸滿足一下條件：A1=16×10-4m2，A2=9×10-4m2，A3=4×10-4m2，Ei=20×1010Pa，i=1,2,3，pi=7.8×103Kg/m3，i=1,2,3，l1=1m，l2=0.5m，l3=0.25m。編寫一個(gè)MATLAB程序解決以下問題。a，在在載荷p3=1000N下u1，u2，u3的位移b，階梯軸的固有頻率和模態(tài)12.8MATLAB應(yīng)用舉例例12.5階梯軸的有限元分析12.8MATLAB應(yīng)用舉例解決方案：階梯軸的剛度矩陣和質(zhì)量矩陣如下所示：12.8MATLAB應(yīng)用舉例解決方案：階梯軸的剛度矩陣和質(zhì)在載荷p3的作用下系統(tǒng)的平衡方程如下所示：