第一章有限單元法的簡要介紹和發(fā)展歷史選編課件

有限單元法

FINITEELEMENTMETHOD

主講：江巍

2015年春季學期有限單元法

FINITEELEMENTMETHOD

參考書目[1]王勖成,邵敏.有限單元法基本原理和數(shù)值方法.清華大學出版社[2]周中堅,盧耀祖.機械與機械結構的有限元分析.同濟大學出版社[3]朱伯芳.有限單元原理及其應用.中國水利水電出版社[4]蔣孝煜.有限元法基礎.清華大學出版社[5]徐芝綸.彈性力學簡明教程.高等教育出版社參考書目[1]王勖成,邵敏.有限單元法基本原理和主要有德國的ASKA；英國的PAFEC；法國的SYSTUS；美國的ALGOR、ABQUS、ADINA、ANSYS、SAP90、BERSAFE、BOSOR、COSMOS、ELAS、MARC和STARDYNE等公司的產(chǎn)品。商業(yè)軟件主要有德國的ASKA；商業(yè)軟件第一章

有限單元法發(fā)展歷史和簡要介紹

第一章

有限單元法發(fā)展歷史和簡要介紹

18世紀末，歐拉在創(chuàng)立變分法的同時就曾用與現(xiàn)代有限元相似的方法求解軸力桿的平衡問題1943年Courant用最小勢能原理和現(xiàn)代有限元法中的線性三角元求解stVenant彈性扭轉問題1952_1853期間，R.W.Clough和M.J.Turner在分析三角翼振動問題時，提出了把平面平面應力三角形板組合起來表達機翼剛度方法，當時稱為直接剛度法。1956年M.J.Turner，R.W.Martin,L.J.Toop在紐約舉行的航空年會上發(fā)表論文《復雜結構的剛度和變形分析》1960年R.W.Clough在論文《平面應力分析的有限單元法》中,首次提出了有限單元,，他因此被稱為“有限單元之父”。《JournalofAppliedMechanics》許多年都拒絕刊登關于有限元方法的文章。發(fā)展歷史之啟蒙18世紀末，歐拉在創(chuàng)立變分法的同時就曾用與現(xiàn)代有限元相似的方

我國已故著名計算數(shù)學專家馮康教授也獨立創(chuàng)立了有限元法，為什么這么說呢？這是由我國當時特定的歷史環(huán)境所決定的。曾經(jīng)有很長一段時間，我國的學術界處于與世隔決的狀態(tài)。正因如此，他的工作才得到了全世界的承認。他最初提出這個方法時，并不知道“有限元”這個名詞，因此他將自己的方法稱之為“基于變分原理的差分格式”。發(fā)展歷史之啟蒙我國已故著名計算數(shù)學專家馮康教授也獨立創(chuàng)立了有眾多數(shù)學家的加盟使得有限元進入黃金發(fā)展階段。有限元方法的理論和程序主要來自各個高校和實驗室Berkeley的EdWilson發(fā)布了第一個程序，第一代的程序沒有名字，第二代線性程序就是著名的SAP(structuralanalysisprogram)，非線性程序就是NONSAP。位于洛杉磯的MSC公司自1963創(chuàng)立并開發(fā)了結構分析軟件SADSAM，在NASA項目資助下MSC于1971年推出自己的專利版本MSC.Nastran。第一批非線性有限元方法的主要貢獻者有Argyris(1965)，Marcal和King(1967)，其中PedroMarcal畢業(yè)于Berkeley大學，任教于Brown大學，于1969年創(chuàng)建了第一家非線性有限元軟件公司MARC公司，在1999年被MSC公司收購。發(fā)展歷史之誕生眾多數(shù)學家的加盟使得有限元進入黃金發(fā)展階段。發(fā)展歷史K.J.Bathe（導師EdWilson），MIT任教，在NONSAP的基礎上發(fā)表了著名的非線性求解器ADINA(AutomaticDynamicIncrementalNonlinearAnalysis)，其源代碼因為長時期廣泛流傳而容易獲得。DavidHibbitt（導師PedroMarcal），在1972年與Karlsson和Sorensen共同建立HKS公司，推出了Abaqus軟件。Abaqus憑借強大的技術、出色的前后處理和可拓展的二次開發(fā)功能，穩(wěn)占高校和研究所的市場，論文發(fā)表數(shù)量多。JohnSwanson博士在Westinghouse公司為核能應用方面發(fā)展了一個非線性有限元程序(主要是關注非線性材料)，于1970年創(chuàng)建SASI(SwansonAnalysisSystem,Inc)公司，后來重組更名為ANSYS公司，ANSYS是著名的多物理材料非線性有限元軟件，通過并購發(fā)展迅速壯大，模塊越來越多，商業(yè)化程度和市場占有率很高。發(fā)展歷史之崛起）K.J.Bathe（導師EdWilson），MIT任教，與其它課程的關系與其它課程的關系各門課程的任務材料力學：研究桿狀構件在拉壓，剪切，彎曲，扭轉作用下的應力和位移。結構力學：在材料力學基礎上研究桿狀構件所組成的結構例如，行架，剛架等，這些都是所謂的桿件系統(tǒng)。彈性力學：非桿狀結構，例如板和水壩，地基等實體結構以及對桿狀構件作進一步，較精確的分析。它與材料力學的研究方法不同，主要是在材力中引入了構件形變狀態(tài)或應力分布的假設，使數(shù)學推導大大簡化，其解是理論解（近似的），而彈性力學則更精確一些。計算力學：是應用結構力學，彈性力學，計算數(shù)學，計算機學的一個結合，提供近似的數(shù)值計算方法，解決問題，而有限元法是其中的一種方法。

上述各種方法最終目標是確立研究對象的應力，形變和位移，用以校核其是否有所需要的強度和剛度。各門課程的任務材料力學：研究桿狀構件在拉壓，剪切，彎曲，扭轉關于有限元法英文縮寫FEM（FiniteElementMethod）應用中習慣稱有限元分析是一種連續(xù)結構離散化數(shù)值計算方法，借助于數(shù)學和力學知識，利用計算機技術而解決工程技術問題FEM與CAECAE－計算機輔助工程（ComputerAidedEngineering）CAE范圍更廣，還包含其它工程分析方法基本思想關于有限元法基本思想基本思想將一個連續(xù)的求解域（連續(xù)體）離散化即分割成彼此用節(jié)點（離散點）互相聯(lián)系的有限個單元，在單元體內(nèi)假設近似解的模式，用有限個結點上的未知參數(shù)表征單元的特性，然后用適當?shù)姆椒?，將各個單元的關系式組合成包含這些未知參數(shù)的代數(shù)方程，得出個結點的未知參數(shù)，再利用插值函數(shù)求出近似解。是一種有限的單元離散某連續(xù)體然后進行求解得一種數(shù)值計算的近似方法。由于單元可以被分割各種形狀和大小不同的尺寸，所以它能很好的適應復雜的幾何形狀，復雜的材料特性和復雜的邊界條件，再加上它有成熟的大型軟件系統(tǒng)支持，使它已成為一種非常受歡迎的，應用極廣的數(shù)值計算方法。基本思想將一個連續(xù)的求解域（連續(xù)體）離散化即分割成彼此用操作流程位移型有限元法求解靜力問題的一般步驟：１）劃分單元；２）計算單元剛度矩陣；３）進行載荷移置；４）引入約束，解方程組求得位移；５）計算應力和應變。注：若以節(jié)點力為未知參數(shù)，先求出節(jié)點處的節(jié)點力，后求位移與應力的方法，稱為力型有限元法。操作流程位移型有限元法求解靜力問題的一般步驟：操作流程結構離散化：

1）劃分網(wǎng)格；

2）載荷移置；

3）簡化約束。單元剛度矩陣與剛度系數(shù)：

1）單元剛度矩陣物理意義為單元抵抗變形的能力；

2）剛度系數(shù)的物理意義是產(chǎn)生單位位移時需要的力的大小。操作流程結構離散化：mm2

，

mm2

，

mmMPa

kN

N試計算應力。

分析過程如下：1．離散化將桿劃分為兩個單元的集合

共有三個節(jié)點

簡單實例mm2，mm2，mmMPakNN試計2．確定單元位移模式（即單元位移函數(shù)）單元e的內(nèi)部，位移按線性規(guī)律變化，即

(1)

本例中每個節(jié)點只有一個自由度，對單元及節(jié)點自由度進行編號簡單實例選擇三個節(jié)點位移

、

、

為基本未知量。2．確定單元位移模式（即單元位移函數(shù)）單元e的內(nèi)部，位移按任取一個單元e作為考察對象，確定位移函數(shù)中系數(shù)a,b在有限元分析過程中，為方便起見，通常使用兩套不同的坐標系。

一是整個結構的參照系oxyz，稱為整體坐標系另一套坐標oxyz

建立在每個單元上，坐標原點和指向都隨單元而變，這種只對單元有效的坐標系，稱為局部坐標系(localcoordinatesystem)。

簡單實例e任取一個單元e作為考察對象，確定位移函數(shù)中系數(shù)a,b在有限元簡單實例簡單實例任意常數(shù)a、b由單元e內(nèi)兩節(jié)點i、j的位移值確定，即：i節(jié)點：

j節(jié)點：

求得：

(2)

代入位移函數(shù)：

為確定系數(shù)a,b，本例使用如圖局部坐標系統(tǒng)簡單實例任意常數(shù)a、b由單元e內(nèi)i節(jié)點：j節(jié)點：求得：得到單元e內(nèi)任意一點x的位移表達式為形狀函數(shù)或形函數(shù)單元節(jié)點位移矢量

簡單實例

(3)

得到單元e內(nèi)任意一點x的位移表達式為形狀函數(shù)或形函數(shù)單元節(jié)點3．推導單元剛度矩陣和單元節(jié)點荷載單元剛度矩陣可由最小勢能原理導出

其中，單元內(nèi)力所做的虛功：整個結構的總勢能

為

簡單實例

(4)

(5)

3．推導單元剛度矩陣和單元節(jié)點荷載單元剛度矩陣可由最小勢能原依據(jù)彈性力學位移與應變的關系得

簡單實例

(6)

依據(jù)彈性力學位移與應變的關系得簡單實例根據(jù)虎克定律：代入單元內(nèi)力虛功表達式，得到將代入上式，得到改寫為簡單實例

(7)

根據(jù)虎克定律：代入單元內(nèi)力虛功表達式，得到將代入上式，得到改積分，得用矩陣寫成

其中，單元剛度矩陣這里

簡單實例

(8)

積分，得用矩陣寫成其中，單元剛度矩陣這里簡單實例單元外力所做虛功為單元節(jié)點荷載簡單實例

(9)

單元外力所做虛功為單元節(jié)點荷載簡單實例在荷載作用下，結構處于平衡狀態(tài)。

則，整個結構的總勢能

為

則由最小勢能原理

，i=1,2,3簡單實例

(10)

在荷載作用下，結構處于平衡狀態(tài)。則，整個結構的總勢能為即：若將最小勢能原理用于單個單元，則得到任一單元的平衡條件為

簡單實例

(11)

即：若將最小勢能原理用于單個單元，簡單實例4．組集總體剛度矩陣和荷載矢量將(8)式中的單元剛度矩陣將式(9)中的單元節(jié)點荷載矩陣組集成整體節(jié)點荷載矩陣最后得到系統(tǒng)的整體平衡方程整體剛度矩陣節(jié)點位移節(jié)點荷載簡單實例組集成整體剛度矩陣

(12)

4．組集總體剛度矩陣和荷載矢量將(8)式中的單元剛度矩陣首先根據(jù)已知數(shù)據(jù)，計算各單元剛度矩陣單元1：簡單實例首先根據(jù)已知數(shù)據(jù)，計算各單元剛度矩陣單元1：簡單實例單元2：簡單實例單元2：簡單實例組集（對號入座）由于總的未知量有三個，所以總剛度矩陣的階數(shù)應為3×3

將中與

對應的元素相加簡單實例組集（對號入座）由于總的未知量有三個，所以總剛度矩陣的階數(shù)應節(jié)點力矢量為：所以，總平衡方程為簡單實例節(jié)點力矢量為：所以，總平衡方程為簡單實例5．約束處理、由于對位移未加任何限制，所以從方程(10.22)得不到節(jié)點位移、的唯一解答。

從數(shù)學上來說，矩陣所以得不到未知的位移分量。

具有奇異性（行列式的值為零），不可求逆，為此，必須引入幾何邊界條件，對方程進行修改。

簡單實例5．約束處理、由于對位移未加任何限制，所以從方程(10.22本例節(jié)點1不能移動，即這里采用

主對角元置1法

即在中，將與零位移的元素改為零，主對角元素改為1；

對應的行和列（本例為第一行、第一列）同時將中與對應的行（未知反力）也改為零。這樣，方程(12)變成：

簡單實例

(13)

本例節(jié)點1不能移動，即這里采用主對角元置1法即在中，將與即

方程(13)的系數(shù)矩陣非奇異，可解。213123簡單實例即方程(13)的系數(shù)矩陣非奇異，可解。213123簡單實例6．求位移解方程(13)，得位移mm

，mm7．求內(nèi)力、應力和應變將(13)式的位移結果代入(11)式，可得單元節(jié)點力矢量簡單實例6．求位移解方程(13)，得位移mm，mm7．求內(nèi)力、工程實例液壓挖掘機有限元分析工程實例工程實例滲流問題有限元分析工程實例工程實例結構問題有限元分析工程實例工程實例隧道問題有限元分析工程實例人有了知識，就會具備各種分析能力，明辨是非的能力。所以我們要勤懇讀書，廣泛閱讀，古人說“書中自有黃金屋。”通過閱讀科技書籍，我們能豐富知識，培養(yǎng)邏輯思維能力；通過閱讀文學作品，我們能提高文學鑒賞水平，培養(yǎng)文學情趣；通過閱讀報刊，我們能增長見識，擴大自己的知識面。有許多書籍還能培養(yǎng)我們的道德情操，給我們巨大的精神力量，鼓舞我們前進。人有了知識，就會具備各種分析能力，第一章有限單元法的簡要介紹和發(fā)展歷史選編課件有限單元法

FINITEELEMENTMETHOD

主講：江巍

2015年春季學期有限單元法

FINITEELEMENTMETHOD

參考書目[1]王勖成,邵敏.有限單元法基本原理和數(shù)值方法.清華大學出版社[2]周中堅,盧耀祖.機械與機械結構的有限元分析.同濟大學出版社[3]朱伯芳.有限單元原理及其應用.中國水利水電出版社[4]蔣孝煜.有限元法基礎.清華大學出版社[5]徐芝綸.彈性力學簡明教程.高等教育出版社參考書目[1]王勖成,邵敏.有限單元法基本原理和主要有德國的ASKA；英國的PAFEC；法國的SYSTUS；美國的ALGOR、ABQUS、ADINA、ANSYS、SAP90、BERSAFE、BOSOR、COSMOS、ELAS、MARC和STARDYNE等公司的產(chǎn)品。商業(yè)軟件主要有德國的ASKA；商業(yè)軟件第一章

有限單元法發(fā)展歷史和簡要介紹

第一章

有限單元法發(fā)展歷史和簡要介紹

18世紀末，歐拉在創(chuàng)立變分法的同時就曾用與現(xiàn)代有限元相似的方法求解軸力桿的平衡問題1943年Courant用最小勢能原理和現(xiàn)代有限元法中的線性三角元求解stVenant彈性扭轉問題1952_1853期間，R.W.Clough和M.J.Turner在分析三角翼振動問題時，提出了把平面平面應力三角形板組合起來表達機翼剛度方法，當時稱為直接剛度法。1956年M.J.Turner，R.W.Martin,L.J.Toop在紐約舉行的航空年會上發(fā)表論文《復雜結構的剛度和變形分析》1960年R.W.Clough在論文《平面應力分析的有限單元法》中,首次提出了有限單元,，他因此被稱為“有限單元之父”?！禞ournalofAppliedMechanics》許多年都拒絕刊登關于有限元方法的文章。發(fā)展歷史之啟蒙18世紀末，歐拉在創(chuàng)立變分法的同時就曾用與現(xiàn)代有限元相似的方

我國已故著名計算數(shù)學專家馮康教授也獨立創(chuàng)立了有限元法，為什么這么說呢？這是由我國當時特定的歷史環(huán)境所決定的。曾經(jīng)有很長一段時間，我國的學術界處于與世隔決的狀態(tài)。正因如此，他的工作才得到了全世界的承認。他最初提出這個方法時，并不知道“有限元”這個名詞，因此他將自己的方法稱之為“基于變分原理的差分格式”。發(fā)展歷史之啟蒙我國已故著名計算數(shù)學專家馮康教授也獨立創(chuàng)立了有眾多數(shù)學家的加盟使得有限元進入黃金發(fā)展階段。有限元方法的理論和程序主要來自各個高校和實驗室Berkeley的EdWilson發(fā)布了第一個程序，第一代的程序沒有名字，第二代線性程序就是著名的SAP(structuralanalysisprogram)，非線性程序就是NONSAP。位于洛杉磯的MSC公司自1963創(chuàng)立并開發(fā)了結構分析軟件SADSAM，在NASA項目資助下MSC于1971年推出自己的專利版本MSC.Nastran。第一批非線性有限元方法的主要貢獻者有Argyris(1965)，Marcal和King(1967)，其中PedroMarcal畢業(yè)于Berkeley大學，任教于Brown大學，于1969年創(chuàng)建了第一家非線性有限元軟件公司MARC公司，在1999年被MSC公司收購。發(fā)展歷史之誕生眾多數(shù)學家的加盟使得有限元進入黃金發(fā)展階段。發(fā)展歷史K.J.Bathe（導師EdWilson），MIT任教，在NONSAP的基礎上發(fā)表了著名的非線性求解器ADINA(AutomaticDynamicIncrementalNonlinearAnalysis)，其源代碼因為長時期廣泛流傳而容易獲得。DavidHibbitt（導師PedroMarcal），在1972年與Karlsson和Sorensen共同建立HKS公司，推出了Abaqus軟件。Abaqus憑借強大的技術、出色的前后處理和可拓展的二次開發(fā)功能，穩(wěn)占高校和研究所的市場，論文發(fā)表數(shù)量多。JohnSwanson博士在Westinghouse公司為核能應用方面發(fā)展了一個非線性有限元程序(主要是關注非線性材料)，于1970年創(chuàng)建SASI(SwansonAnalysisSystem,Inc)公司，后來重組更名為ANSYS公司，ANSYS是著名的多物理材料非線性有限元軟件，通過并購發(fā)展迅速壯大，模塊越來越多，商業(yè)化程度和市場占有率很高。發(fā)展歷史之崛起）K.J.Bathe（導師EdWilson），MIT任教，與其它課程的關系與其它課程的關系各門課程的任務材料力學：研究桿狀構件在拉壓，剪切，彎曲，扭轉作用下的應力和位移。結構力學：在材料力學基礎上研究桿狀構件所組成的結構例如，行架，剛架等，這些都是所謂的桿件系統(tǒng)。彈性力學：非桿狀結構，例如板和水壩，地基等實體結構以及對桿狀構件作進一步，較精確的分析。它與材料力學的研究方法不同，主要是在材力中引入了構件形變狀態(tài)或應力分布的假設，使數(shù)學推導大大簡化，其解是理論解（近似的），而彈性力學則更精確一些。計算力學：是應用結構力學，彈性力學，計算數(shù)學，計算機學的一個結合，提供近似的數(shù)值計算方法，解決問題，而有限元法是其中的一種方法。

上述各種方法最終目標是確立研究對象的應力，形變和位移，用以校核其是否有所需要的強度和剛度。各門課程的任務材料力學：研究桿狀構件在拉壓，剪切，彎曲，扭轉關于有限元法英文縮寫FEM（FiniteElementMethod）應用中習慣稱有限元分析是一種連續(xù)結構離散化數(shù)值計算方法，借助于數(shù)學和力學知識，利用計算機技術而解決工程技術問題FEM與CAECAE－計算機輔助工程（ComputerAidedEngineering）CAE范圍更廣，還包含其它工程分析方法基本思想關于有限元法基本思想基本思想將一個連續(xù)的求解域（連續(xù)體）離散化即分割成彼此用節(jié)點（離散點）互相聯(lián)系的有限個單元，在單元體內(nèi)假設近似解的模式，用有限個結點上的未知參數(shù)表征單元的特性，然后用適當?shù)姆椒?，將各個單元的關系式組合成包含這些未知參數(shù)的代數(shù)方程，得出個結點的未知參數(shù)，再利用插值函數(shù)求出近似解。是一種有限的單元離散某連續(xù)體然后進行求解得一種數(shù)值計算的近似方法。由于單元可以被分割各種形狀和大小不同的尺寸，所以它能很好的適應復雜的幾何形狀，復雜的材料特性和復雜的邊界條件，再加上它有成熟的大型軟件系統(tǒng)支持，使它已成為一種非常受歡迎的，應用極廣的數(shù)值計算方法?；舅枷雽⒁粋€連續(xù)的求解域（連續(xù)體）離散化即分割成彼此用操作流程位移型有限元法求解靜力問題的一般步驟：１）劃分單元；２）計算單元剛度矩陣；３）進行載荷移置；４）引入約束，解方程組求得位移；５）計算應力和應變。注：若以節(jié)點力為未知參數(shù)，先求出節(jié)點處的節(jié)點力，后求位移與應力的方法，稱為力型有限元法。操作流程位移型有限元法求解靜力問題的一般步驟：操作流程結構離散化：

1）劃分網(wǎng)格；

2）載荷移置；

3）簡化約束。單元剛度矩陣與剛度系數(shù)：

1）單元剛度矩陣物理意義為單元抵抗變形的能力；

2）剛度系數(shù)的物理意義是產(chǎn)生單位位移時需要的力的大小。操作流程結構離散化：mm2

，

mm2

，

mmMPa

kN

N試計算應力。

分析過程如下：1．離散化將桿劃分為兩個單元的集合

共有三個節(jié)點

簡單實例mm2，mm2，mmMPakNN試計2．確定單元位移模式（即單元位移函數(shù)）單元e的內(nèi)部，位移按線性規(guī)律變化，即

(1)

本例中每個節(jié)點只有一個自由度，對單元及節(jié)點自由度進行編號簡單實例選擇三個節(jié)點位移

、

、

為基本未知量。2．確定單元位移模式（即單元位移函數(shù)）單元e的內(nèi)部，位移按任取一個單元e作為考察對象，確定位移函數(shù)中系數(shù)a,b在有限元分析過程中，為方便起見，通常使用兩套不同的坐標系。

一是整個結構的參照系oxyz，稱為整體坐標系另一套坐標oxyz

建立在每個單元上，坐標原點和指向都隨單元而變，這種只對單元有效的坐標系，稱為局部坐標系(localcoordinatesystem)。

簡單實例e任取一個單元e作為考察對象，確定位移函數(shù)中系數(shù)a,b在有限元簡單實例簡單實例任意常數(shù)a、b由單元e內(nèi)兩節(jié)點i、j的位移值確定，即：i節(jié)點：

j節(jié)點：

求得：

(2)

代入位移函數(shù)：

為確定系數(shù)a,b，本例使用如圖局部坐標系統(tǒng)簡單實例任意常數(shù)a、b由單元e內(nèi)i節(jié)點：j節(jié)點：求得：得到單元e內(nèi)任意一點x的位移表達式為形狀函數(shù)或形函數(shù)單元節(jié)點位移矢量

簡單實例

(3)

得到單元e內(nèi)任意一點x的位移表達式為形狀函數(shù)或形函數(shù)單元節(jié)點3．推導單元剛度矩陣和單元節(jié)點荷載單元剛度矩陣可由最小勢能原理導出

其中，單元內(nèi)力所做的虛功：整個結構的總勢能

為

簡單實例

(4)

(5)

3．推導單元剛度矩陣和單元節(jié)點荷載單元剛度矩陣可由最小勢能原依據(jù)彈性力學位移與應變的關系得

簡單實例

(6)

依據(jù)彈性力學位移與應變的關系得簡單實例根據(jù)虎克定律：代入單元內(nèi)力虛功表達式，得到將代入上式，得到改寫為簡單實例

(7)

根據(jù)虎克定律：代入單元內(nèi)力虛功表達式，得到將代入上式，得到改積分，得用矩陣寫成

其中，單元剛度矩陣這里

簡單實例

(8)

積分，得用矩陣寫成其中，單元剛度矩陣這里簡單實例單元外力所做虛功為單元節(jié)點荷載簡單實例

(9)

單元外力所做虛功為單元節(jié)點荷載簡單實例在荷載作用下，結構處于平衡狀態(tài)。

則，整個結構的總勢能

為

則由最小勢能原理

，i=1,2,3簡單實例

(10)

在荷載作用下，結構處于平衡狀態(tài)。則，整個結構的總勢能為即：若將最小勢能原理用于單個單元，則得到任一單元的平衡條件為

簡單實例

(11)

即：若將最小勢能原理用于單個單元，簡單實例4．組集總體剛度矩陣和荷載矢量將(8)式