有限元分析程序設(shè)計(jì)

1、精選優(yōu)質(zhì)文檔-傾情為你奉上 結(jié)構(gòu)有限元分析程序設(shè)計(jì)緒 論§0.1 開設(shè)“有限元程序設(shè)計(jì)”課程的意義和目的§0.2 課程特點(diǎn)§0.3 課程安排§0.4 課程要求§0.5 基本方法復(fù)習(xí)$0.1 意義和目的1. 有限元數(shù)值分析技術(shù)本身要求工程設(shè)計(jì)研究人員掌握1). 有限元數(shù)值分析技術(shù)的完善標(biāo)志著現(xiàn)代計(jì)算力學(xué)的真正成熟和實(shí)用化,已在各種力學(xué)中得到了廣泛的應(yīng)用。比如： ，已楊為工程結(jié)構(gòu)分析中最得以收斂的技術(shù)手段，現(xiàn)代功用大致有：a) . 現(xiàn)代結(jié)構(gòu)論證。對結(jié)構(gòu)設(shè)計(jì)從內(nèi)力，位移等方面進(jìn)行優(yōu)劣評定，從而進(jìn)行結(jié)構(gòu)優(yōu)化設(shè)計(jì)。b) 可取代部份實(shí)驗(yàn)，局部實(shí)驗(yàn)+有限元分

2、析，是現(xiàn)代工程設(shè)計(jì)研究方法的一大特點(diǎn)。c) 結(jié)構(gòu)的各種功能分析（疲勞斷裂，可靠性分析等）都以有限元分析工具作為核心的計(jì)算工具。 2). 有限元數(shù)值分析本身包括著理論+技術(shù)實(shí)現(xiàn)（本身功用所絕定的） 有限元數(shù)值分析本身包括著泛函理論+分片插值函數(shù)+程序設(shè)計(jì)2. 有限元分析的技術(shù)實(shí)現(xiàn)（近十佘年的事）更依賴于計(jì)算機(jī)程序設(shè)計(jì)有限元分析的技術(shù)取得的巨大的成就，從某種意義上說，得益于計(jì)算機(jī)硬件技術(shù)的發(fā) 展和程序設(shè)計(jì)技術(shù)的發(fā)展，這兩者的依賴性在當(dāng)代表現(xiàn)得更加突出。（如可視化技術(shù)）3. 從學(xué)習(xí)的角度，不僅要學(xué)習(xí)理論，而且要從程序設(shè)計(jì)設(shè)計(jì)角度對這些理論的技術(shù)實(shí)現(xiàn)有一個深入的了解，應(yīng)當(dāng)致力于掌握這些技術(shù)實(shí)現(xiàn)能力，從

3、而開發(fā)它，發(fā)展它。（理論本身還有待于進(jìn)一步完美相應(yīng)的程序設(shè)計(jì)必須去開發(fā)）4. 程序設(shè)計(jì)不僅是實(shí)現(xiàn)有限元數(shù)值分析的工具和橋梁，而且在以下諸方面也有意義：1). 精通基本概念，深化理論認(rèn)識；2). 鍛煉實(shí)際工程分析，實(shí)際動手的能力；3). 獲得以后工作中必備的工具。（作業(yè)+老師給元素庫）目的：通過講述有限元程序設(shè)計(jì)的技術(shù)與技巧，便能達(dá)到自編自讀的能力。§0.2 課程特點(diǎn) 總描述： 理論+算法+數(shù)據(jù)結(jié)構(gòu)（程序設(shè)計(jì)的意義） 理論： 有限元算法，構(gòu)造，步驟，解的等外性，收斂性，穩(wěn)定性，誤差分析 算法； 指求解過程的技術(shù)方法，含兩方面的含義；a. 有限元數(shù)值分析算法，b, 與數(shù)據(jù)結(jié)構(gòu)有關(guān)的算法（

4、總剛稀疏存貯，提取，節(jié)點(diǎn)優(yōu)化編號等） 數(shù)據(jù)結(jié)構(gòu)： 指各向量矩陣存貯管理與實(shí)現(xiàn)，輔助管理結(jié)構(gòu)（指針，數(shù)據(jù)記錄等） 具體特點(diǎn)：理論性強(qiáng)： 能量泛函理論+有限元構(gòu)造算法+數(shù)據(jù)結(jié)構(gòu)構(gòu)造算法內(nèi)容繁雜： 理論方法+技術(shù)方法+技術(shù)技巧技巧性強(qiáng)： 排序，管理結(jié)構(gòu)（指針生成，整型運(yùn)算等）§0.3 課程安排. 單元剛度矩陣及元素設(shè)計(jì)（單元剛陣算法，桿梁平面分析，板彎非協(xié)調(diào)元等）. 總剛的形式及程序設(shè)計(jì)（單剛提前準(zhǔn)備，技術(shù)復(fù)雜）. l邊界條件及程序設(shè)計(jì)（等效荷載計(jì)算，位移邊界條件置入，多工況的對稱性）. 總剛線性方程組求解（LDLT分解，分塊算法，子結(jié)構(gòu)算法，波前法） 單元應(yīng)力計(jì)算+應(yīng)力處理與改善。. 數(shù)

5、據(jù)處理（數(shù)據(jù)分類，壓縮存貯，節(jié)點(diǎn)排序方法）, 變帶寬存貯的優(yōu)化理論，圖論的理論，有限元的圖結(jié)構(gòu)，存貯管理復(fù)核）. 有限元議程全稀疏管理與求解策略。說明： 僅線性部份，復(fù)材，接觸，彈塑性等不包括，基本部份。實(shí)踐性作業(yè)安排：1. 作業(yè)： 總的結(jié)構(gòu)管理程序+子功能模塊的編程，一個題的計(jì)算實(shí)踐2. 送有限元元素庫。§0.4 課程要求1. 先行要求2. 作業(yè)要求（計(jì)算機(jī)編程得出正確結(jié)果）3. 課程校核要求§0.5 基本方法復(fù)習(xí)0.5.1 結(jié)構(gòu)化程序設(shè)計(jì)方法0.5.2 有限元分析方法回顧0.5.3 Fortran語言回顧0.5.4 結(jié)構(gòu)化程序設(shè)計(jì)1. 基本結(jié)構(gòu)：構(gòu)成一個問題從輸入到求解

6、輸出的基本程序形式：ImputProcessOutput(輸入) （處理） （輸出）三種基本形式：a). 樹形結(jié)構(gòu)（順序執(zhí)行結(jié)構(gòu)） Imputp1 . pn Output (多用于程序結(jié)構(gòu) call1,call2.)b). t選擇結(jié)構(gòu)（條件執(zhí)行結(jié)構(gòu)） 內(nèi)部算法實(shí)現(xiàn): IF.THEN; SWHICH,CASEc).循環(huán)結(jié)構(gòu)(反復(fù)執(zhí)行結(jié)構(gòu))特點(diǎn): 結(jié)構(gòu)特征簡單明了,易讀,易調(diào)試. 盡量少用GO TO 語句2. 整體結(jié)構(gòu)(算法語言系統(tǒng)結(jié)構(gòu))積木式(Fortran): 每個設(shè)定的功能分析團(tuán)體的一個模塊,每個模塊又稱作整體結(jié)構(gòu)的素材,主模塊象積木一樣堆積.語體不聯(lián)系,但有通訊方法溝通模塊間信息,各種模塊

7、有各自的特征語體,main progame,subroutine sub. 嵌套式(Pascall): 主模塊與子模塊相互嵌套,各模塊的特征,語體相同 procedure main procefure ssub1 . procedure sub2 . End sub1 end sub2 . End main 函數(shù)式(c語言): 主要特點(diǎn)是功能模塊作為庫函數(shù)調(diào)用,需用時在庫內(nèi)調(diào)用,每一個函數(shù)有表征語句,這種語體接近自然思維,而且對系統(tǒng)資源的調(diào)配應(yīng)用更完善. 面向?qū)ο蟮某绦蛟O(shè)計(jì): 實(shí)施過程的可視化+控制性 3. 結(jié)構(gòu)化程序設(shè)計(jì)方法a). TopDown(自上面下),系統(tǒng)性強(qiáng),選擇性強(qiáng).b). Cr

8、itical Component First (關(guān)鍵部份優(yōu)先),先抓主要矛盾,分清重,緩,急.c). 獨(dú)立調(diào)試,總體聯(lián)調(diào),(軟件設(shè)計(jì)的社會化作業(yè)).4. 程序設(shè)計(jì)要點(diǎn) a). 自覺有意識地設(shè)計(jì)一個良好的程序結(jié)構(gòu),做到:易讀,易懂,易管理,易修改,易發(fā)展. b). 做到邏輯清晰,說明完整. c). 要有工藝設(shè)計(jì)概念有框圖,有步驟.5). 結(jié)構(gòu)化程序設(shè)計(jì)原則 a). 盡可能通用性好（適應(yīng)各種規(guī)模的復(fù)題，？的擴(kuò)大依據(jù)程序設(shè)計(jì)指標(biāo)而定） b).整體精煉，清晰；避免GOTO。c).省機(jī)時，省存貯，計(jì)算精度高，（算法上下功夫，要理論分析加技巧）d).輸入數(shù)據(jù)少，格式簡單。e).輸出結(jié)果簡明，忌諱打印過多（

9、與具體調(diào)試過程不一樣）。f).易讀易維護(hù)，易發(fā)展。§0.5.2.有限無方法求解過程回顧一. 力學(xué)模型的分級管理有限無程序?qū)αW(xué)模型的數(shù)據(jù)按一級：結(jié)構(gòu)級（有點(diǎn)廣義，不僅指具體結(jié)構(gòu)，也指模型題目的規(guī)模） 二級：單元級三級：節(jié)點(diǎn)級2.基本關(guān)系3. 描述參數(shù) A). 節(jié)點(diǎn)描述參數(shù) (1). 節(jié)點(diǎn)位置(總體坐標(biāo)系下的坐標(biāo)). (2). 節(jié)點(diǎn)局部坐標(biāo)(按節(jié)點(diǎn)的約束方向制定的特殊坐標(biāo)系x，y，z,v如斜支撐) (3). 節(jié)點(diǎn)的性質(zhì)(自由,固定,指定位移,從屬其它節(jié)點(diǎn)). (4). 節(jié)點(diǎn)力:(Fx ,Fy , Fz ,Mx ,My ,Mz ) (5). 節(jié)點(diǎn)位移:(u , v , w, x ,y ,

10、z )B. 單元描述參數(shù) (1). 材料特性參數(shù)不清 E, G, D (2). 節(jié)點(diǎn)的幾何剛度參數(shù)(即面積A,板厚H,梁抗彎模量I) (3). 單元的局部坐標(biāo). (用于應(yīng)力分析等,如圖形曲面) (4). 單元的節(jié)點(diǎn)編號 (5). 單元的幾何矩陣營 (節(jié)點(diǎn)變形與應(yīng)力關(guān)系矩陣) (6). 單元剛度矩陣 K (7). 單元的應(yīng)力,應(yīng)變向量,(有限元分析多用向量,而不用矩陣(張量)結(jié)構(gòu)描述參數(shù) 單元總數(shù),節(jié)點(diǎn)總數(shù),單元婁型總數(shù),結(jié)構(gòu)材料種婁數(shù),節(jié)點(diǎn)自由度數(shù)(控制題目規(guī)模)二. 基本公式系統(tǒng)1. 單元剛度計(jì)算公式 2. 單元剛陣組合 K=ATKA 3. 單元節(jié)點(diǎn)荷載計(jì)算4. 節(jié)點(diǎn)荷載組裝: 5. 位移約

11、束關(guān)系:6. 總剛方程解:7. 應(yīng)變計(jì)算:8. 應(yīng)力計(jì)算:9. 支撐反力計(jì)算:輸入邊界條件(對稱條件)形成各荷載工況的節(jié)點(diǎn)荷載陣 總剛分解 回代求出位移及輸出 計(jì)算應(yīng)變、應(yīng)力形成單元剛陣單剛向總剛投放坐標(biāo)變換輸入原始參數(shù)計(jì)算總剛規(guī)模形成總剛方程向總節(jié)點(diǎn)荷載陣投放形成單元荷載陣調(diào)整幾何、彈性矩陣調(diào)整單元位移列陣三. 有限元分析的模塊組織.四 結(jié)構(gòu)分析的原始輸入數(shù)據(jù) 1題目規(guī)模 節(jié)點(diǎn)數(shù)目：NNP 單元數(shù)目： NE 2節(jié)點(diǎn)數(shù)據(jù) 單元人坐標(biāo)：XE（NNP，3） 3單元數(shù)據(jù) 單元節(jié)點(diǎn)編號 ： ME（NE，3）、ME（NE，2） 材料特性： E、N 單元幾何參數(shù)： I、RI（慣矩） 4荷載數(shù)據(jù) 外荷載作用

12、點(diǎn)，坐標(biāo)及大小： PA（NNP，1）§0. 5.3 Fortran語言回顧1子模塊（子程序） subroutinea. 特點(diǎn)：獨(dú)立性強(qiáng)，只要輸入輸出接口，象一個黑匣子，與外界無關(guān)。b. 作用：完成一個獨(dú)立的功能（求應(yīng)力，矩陣分解，投放等）c. 格式： subroutine function(ip1,ip2,rp1,rp2,io1,io2,ro1,ro2).(其中ip1,ip2,rp1,rp2,是輸入形參，io1,io2.rp1,rp2是輸出形參)2數(shù)據(jù)傳遞形式 1）. COMMON 公共塊語句傳遞，（公共塊的內(nèi)容不能作為形參）a. 公共塊分為無名公共塊和有名公共塊b. 公共塊的參數(shù)不

13、能作為子程序的參數(shù)出現(xiàn)，c. 公共塊名一致，其內(nèi)容在不同公共塊中可以標(biāo)志符不同（但其長度應(yīng)一致）d. 通用原始數(shù)據(jù)放入公共塊（作為實(shí)參錯誤率大）e. 盡可能不放數(shù)值，安息組一般可作成可調(diào)長度 f. 格式 Common/comm/. Subroutine fun() Common/comm/.2).形參實(shí)參對應(yīng)a. 實(shí)參不能開辟存貯單元，子程序內(nèi)定義語句中的形參數(shù)組由主程序定義，在子程序中僅形式定義（即僅說明是數(shù)組，因而大小無所謂）b. 格式：Dimension RP(1000),RO(1000) . Call sub1(RP,RO) END DIMENTION IBANK SUBROUTINE

14、 SUB1(RP,RO,NE)DIMENTION RP(1),RO(1),SP(50)DIMENTION RP(NE,1),RO(NE,1) (形參的動態(tài)定義，實(shí)參不能)3）. 數(shù)組長度自動調(diào)整方法。 PROGRAM MAIN INPLICIT REAL*8 (AH,O-Z) CHARACTAR*20 TR COMMON/COMM/. DIMENTION IBANK( ),RBANK( ),IP1( ),IP2( ) IP1(1)= IP(N)=. IP2(1)= IP2(N)=. CALL SUB1(IBANK(IP1(1),IBANK(IP1(N),RBANK(IP2(1),.) . EN

15、D SUBROUTINE SUB1(II1,IO2,.RI1,.RO1.NE) DIMENTION II1(1) ,IO2(NE,1), RI1(1), RO1(1)*.公共塊分無名與有名公共塊：公共塊參數(shù)不能子程序形參中出現(xiàn)。公共塊一致，其內(nèi)容在不同程序中可以標(biāo)識符不同。（但長度要一致）*.公共塊使用原則：1) 通用原始數(shù)據(jù)代入公共塊，作為實(shí)參出錯率大2) 盡可能不放數(shù)組，數(shù)組要作成可調(diào)的（自動可調(diào)）。數(shù)組長度自動調(diào)整方法。§1.單元剛陣幾元素設(shè)計(jì)主講內(nèi)容：工程常用元素的單元剛陣計(jì)算與編程。§1.1桿、梁的單元剛陣§1.2平面問題的三角形、矩形單元剛陣§

16、;1.3板彎問題的三角形、矩形單元剛陣§1.4非、擬協(xié)調(diào)之介紹及分片檢驗(yàn)方法§1.0本講內(nèi)容序§1.0.1單元剛陣在總剛中的地位和作用由有限元分析過程可以知道，最主要的兩步驟：1) 進(jìn)行單元的特性分析，建立單元級的剛度方程（平衡方程），得到單元的剛度矩陣；注：中各元素的位置與和要有嚴(yán)格的對應(yīng)關(guān)系2) 進(jìn)行結(jié)構(gòu)整體分析，集合所有單元的剛度方程建立全結(jié)構(gòu)的剛度方程，設(shè)得到全結(jié)構(gòu)的剛度矩陣。即A為布爾矩陣總剛由單元剛陣按節(jié)點(diǎn)編號順序投放而成，故單元剛陣是整個有限元方法的第一步；是總裝剛陣的原材料。解釋：單元節(jié)點(diǎn)平衡；結(jié)構(gòu)總剛對號和迭加的原理§1.0.2單元剛陣

17、形式1. 基本公式：（1.1）成為幾何矩陣，作用是把應(yīng)變與位移聯(lián)系起來反映了單元的形變特性與性質(zhì)。固而是坐標(biāo)的函數(shù)，規(guī)模與節(jié)點(diǎn)數(shù)目（自由度數(shù)目）有關(guān)，不同類型的單元有不同的。（依據(jù)設(shè)定的位移形函形式與節(jié)點(diǎn)位移的性質(zhì)）。對常應(yīng)元是常數(shù)。為彈性矩陣，只取決于材料特性E，與坐標(biāo)無關(guān)，對于平面應(yīng)力問題： （1.2）對于平面應(yīng)變問題，將上式Eàà代入上式即得。2. 形成方法1) 直接賦值。有了積分出來的 各元素顯式，即可逐項(xiàng)（元素）直接賦值。在賦值時，有數(shù)字變化規(guī)律的應(yīng)按特點(diǎn)賦值，以縮短程序長度，提高執(zhí)行效率。2) 矩陣乘法。用矩陣乘法進(jìn)行賦值。§1.1桿、梁的單元剛陣&#

18、167;1.1.1桿元xy對于各種變截面等軸力桿其單元剛陣顯示：uj(Pj)y（1.3）為等效截面面積。ui(Pi)x對于不同的截面變化情況可由材料力學(xué)方法計(jì)算得到。*對于等截面桿，可設(shè)定形狀函數(shù)，利用公式（1.1）計(jì)算，即由：,(1.4) (1.5)*對于變截面桿，采用假設(shè)形狀函數(shù)的方法不妥（誤差大），因此，常用材料力學(xué)推倒方法：即公式：à (等軸力)（1.6）相對于1點(diǎn)：同理，相對于2點(diǎn)：于是得：（1.6）令：故得：（1.3）式思考題：變截不均的彎軸力桿的剛陣形式如何？對于滿足桿元截面積隨坐標(biāo)單調(diào)變化的函數(shù)，Aeq計(jì)算見P2830*形成方法可采用直接賦值語句完成見P30FAKE1

19、§1.1.2梁元應(yīng)用公式（1.1）進(jìn)行直接積分，得到顯式，然后用上賦值語句編程。1. 梁元形 及剛陣顯式工程上使用的梁元一般為復(fù)合受力狀態(tài)，包括軸向、扭轉(zhuǎn)及各平面內(nèi)的y彎曲，因此對于一般的三維空間中的梁元節(jié)點(diǎn)位移包括：線性位移、及轉(zhuǎn)角位移、；節(jié)點(diǎn)力 位包括：三個方向的力、和三個方向的彎矩。yxzxz故梁元每個節(jié)點(diǎn)有6個自由度，即：à1) 對于軸向拉伸與扭轉(zhuǎn)可采用線性形函描述即：(1.7)應(yīng)變定義（1.8）截面極慣矩（1.9）注意：這里的體積分、面積分為常數(shù)包括在彈性陣內(nèi)(1.10)2) 在xy平面內(nèi)的彎曲，不計(jì)剪切，節(jié)點(diǎn)位移為；節(jié)點(diǎn)力為形 的選取是依據(jù)彎曲理論（材料力學(xué)）的

20、基本假設(shè)的兩直法線公設(shè)所確定即：（1.11）進(jìn)行積分，可得三次撓度曲線（1.12）代入節(jié)點(diǎn)位移系統(tǒng)可得到形狀函數(shù)（1.13）（1.14）*應(yīng)變定義：梁的曲率定義為廣義定義應(yīng)變，梁的彎矩定義為廣義應(yīng)力，即：（1.15）按梁的撓度方程，可有廣義的應(yīng)力應(yīng)變關(guān)系(1.16)其中代入（1.13式），得：（1.17）（1.18）由式（1.1），于是：（1.19）（1.20）3) 在yz平面內(nèi)的彎曲，節(jié)點(diǎn)位移與力系統(tǒng)為：形函推導(dǎo)過程類似。剛陣顯式：（1.21）4) 綜合上述三種單獨(dú)剛度，可以得空間梁元素在局部坐標(biāo)系下復(fù)合狀態(tài)的剛陣顯示見P4243間，即書公式（224），（注意節(jié)點(diǎn)位移、節(jié)點(diǎn)力的排序）。2.

21、形成梁元剛陣的子程序設(shè)計(jì)方法：直接賦值。見P4344FAKE4§1.2平面問題的三角元、矩形單元剛陣§1.2.1平面應(yīng)力三角元。（常應(yīng)變元）1. 公式。(1.22)Qy(1.23)m其中：在板平面內(nèi)TjR轉(zhuǎn)換123iPx2. 形成B陣的輪換循環(huán)賦值方法。I=1,3I1=I+1-3*(I/3)I2=I+2-3*(I+1)/3)KI=ME(K,I1)KJ=ME(K,I2)B(1,2*I-1)=XE(KI,2)-XE(KJ,2)B的第1行1,3,5之元素B(2,2*I)=XE(KJ,1)-XE(KI,1)B的第2行2,4,6之元素B(3,2*I-1)=B(2,2*I)PP34B(

22、3,2*I)=B(1,2*I-1)未完當(dāng)平面元處于空間位置時，（即節(jié)點(diǎn)坐標(biāo)為三維數(shù)據(jù)）。平面內(nèi)的幾何參數(shù)需要計(jì)算，計(jì)算公式：y1） PQ長度：mQT2） PQ的方向余弦：RjiPx3） PT長度（按線段投影公式）4） TR長度：5） QT長度：6） TR的方向余弦7） 面積8） 剛度陣中的參數(shù)9） 在中各長度的取值記號規(guī)則：利用上面的子程序，可直接賦值形成局部坐標(biāo)下的陣。3. 局部系下剛陣計(jì)算與編程。4（1.24）31mP 其中C為板厚利用矩陣乘法可直接運(yùn)算得到間見P35§1.2.2平面應(yīng)力四節(jié)點(diǎn)單元剛度i1 單元平面內(nèi)局部下形狀函數(shù)。2j節(jié)點(diǎn)位移和節(jié)點(diǎn)力系統(tǒng)：形狀函數(shù)：（1.25）

23、（1.26）（1.27）顯見得不再是常數(shù)陣，說明應(yīng)變在平面內(nèi)線性變化。2.矩陣平面應(yīng)力剛陣形成公式：由于不再是 常數(shù)陣，只能手工積分每一元素，設(shè)得到剛陣顯式，然后用賦值方法編程。見公式PP38-39（2-19）剛度形成子程序FAKE3見PP39-40思考題：變軸力、變截面桿單元剛陣的形成方法。（用柔度方法）作業(yè)：（1）閱讀有限元PP27-60內(nèi)容（2）有一空間梁元，截面尺寸：長度1m. 編程計(jì)算單元剛度矩陣。§1.3.1板彎問題的三角元、矩形單元剛陣作為作業(yè)、閱讀理解，不作掌握要求 § 1.4 協(xié)調(diào)、非協(xié)調(diào)、廣義協(xié)調(diào)及分片檢驗(yàn)1、4、0 引 以有限元數(shù)值分析的技術(shù)實(shí)現(xiàn)為目的

24、本門課程，不僅要求學(xué)生能夠進(jìn)行實(shí)際的工程運(yùn)算；另一方面也需要對解的收斂及精確性有所了解，是能從細(xì)節(jié)計(jì)算到理論性質(zhì)都有所把握，這樣，才能做到全面深入有助于對解結(jié)果得理論分析，此為基本之目的。1、4、1 協(xié)調(diào)、非協(xié)調(diào)介紹 位移法有限元以Ritz的結(jié)構(gòu)最小有限元為基礎(chǔ)，該原理在數(shù)學(xué)上是一個泛函極值（變分）問題，系統(tǒng)勢能可以表為以下數(shù)學(xué)形式：面力勢能內(nèi)力勢能體力勢能=1/2 - - (1) 0 。表述為：在所有滿足內(nèi)部連續(xù)性和運(yùn)動學(xué)邊界條件的位移中滿足平衡方程的位移使系統(tǒng)勢能取駐值。如果駐值是極小點(diǎn)的，則平分行是穩(wěn)定階。又：對于精確于問題的位移函數(shù)，系統(tǒng)勢能的變分可求得關(guān)于問題應(yīng)滿足的所有微分方程：平

25、衡方程邊界條件（幾何關(guān)系及物理方程是自然滿足的）遺憾的是精確位移難得尋找，故一般采用泛函的極小化序列逼近方法。類似于傅立葉級數(shù)逼近函數(shù)那樣，把無窮維空間用有限空間去逼近。在有限元當(dāng)中，當(dāng)元素尺寸趨近于0時（即節(jié)點(diǎn)數(shù)目或節(jié)點(diǎn)自由度數(shù)趨于時），最后的解答若能無限逼近準(zhǔn)確解，那么這樣的位移函數(shù)（或形狀函數(shù)）就稱為收斂的，因此從收斂性及算收斂速度方面提出幾點(diǎn)對形狀函數(shù)的要求：、函數(shù)本身及其導(dǎo)數(shù)應(yīng)在元素上連續(xù)，并含有常數(shù)部分；、元素之間的位移協(xié)調(diào)，不僅節(jié)點(diǎn)處的位移應(yīng)當(dāng)協(xié)調(diào)，沿整個內(nèi)邊界上的位移也應(yīng)當(dāng)協(xié)調(diào)（或稱相容 ）。 、多項(xiàng)式的項(xiàng)數(shù)越多越好，因用高次比低次多項(xiàng)式收斂快。、含有剛體位移（平動包含常數(shù)項(xiàng)，

26、轉(zhuǎn)動包含線性項(xiàng)）。協(xié)調(diào)之： 即滿足、條件的形狀函數(shù)的元素，當(dāng)然能滿足3) 4)條件協(xié)調(diào) 元的收斂率就更高。 協(xié)調(diào)元的性質(zhì)：1) 能夠以單調(diào)趨勢逼近于正確解。如曲線.2) 勢能總是大于最小狀態(tài)，故解得上界。3) 近似剛度k偏大，即元素偏“硬”。4) 近似的位移偏小，即求得位移的下界。能夠以單調(diào)趨勢逼近于正確解。如曲線.勢能總是大于最小狀態(tài)，故解得上界。給點(diǎn)數(shù)近似剛度k偏大，即元素偏“硬”。近似的位移偏小，即求得位移的下界非協(xié)調(diào)元：在彈性力學(xué)中，如板彎曲，相鄰元素不僅要求位移本身連續(xù)，而且要求位移的導(dǎo)數(shù)連續(xù)（板彎邊界上的相容性）。而在工程上能夠保證導(dǎo)數(shù)相容的形變往往難以找到，以致工程上只能采用違反

27、相容原則的一些形狀函數(shù)，由違反相容原則的形函所構(gòu)成的元素稱為非協(xié)調(diào)元。非協(xié)調(diào)元性質(zhì)：、不能以最小位能原理作為它的理論基礎(chǔ)。、解的趨勢可能收斂，可能不收斂，（取決于網(wǎng)格劃分）。對于收斂的趨勢也未必滿足單調(diào)性?？赡苁諗壳€如圖中。QQmjik例： 三角板元： 節(jié)點(diǎn)位移系統(tǒng)： 位移形函：W=CL+CLL+CLL+CL+CLL+CLL+CL+ CLL+CLL+CLLL L =( a+ bx + dy )/2a= xy - xy ; b= y- y ; d= - x+ x 。代入節(jié)點(diǎn)位移參數(shù)：可得：w= N N 變換一下寫法：w= H HH H H H H HH H H = L(3 2L) 7LLL H

28、 = L ( dL dL ) + ( d d )LLL H = L ( bL bL ) + ( b b) LLL H = 27 LLLi +2=k i j k 下標(biāo)按循環(huán)計(jì)算 上述元素只能在結(jié)構(gòu)上做到位移導(dǎo)數(shù)連續(xù)，在邊界上其他點(diǎn)處，位移的法向?qū)?shù)并不連續(xù)，這因?yàn)椋河捎诜ㄏ驅(qū)?shù)是一個完全的二次多項(xiàng)式，在元素的每條邊上，其變化規(guī)律位一條二次拋物線，需要三個點(diǎn)上法向?qū)?shù)的相等條件才能維一確定，故相鄰兩條曲線一般不全重合。故所舉三角板彎元為非協(xié)調(diào)元。例 書P的矩形元，由于坐標(biāo)的交叉雙乘積（不完備），可發(fā)現(xiàn)不該是w 或其導(dǎo)數(shù)， 都是連續(xù)的，這樣只要節(jié)點(diǎn)的這些參數(shù)相同，邊界上的這些是沒有問題的，但展開&#

29、160;N的項(xiàng)，可以發(fā)現(xiàn)xy項(xiàng)，或者說缺少了代表熱率變形的一項(xiàng)，因此，作為形狀函數(shù)，是不能保證向正確的解答收斂，因而是非協(xié)調(diào)元。改進(jìn)方案之一，是在節(jié)點(diǎn)處增加節(jié)點(diǎn)參數(shù)，并采用完全的埃爾米特三次多項(xiàng)式。§1、4、2 非協(xié)調(diào)元的排先檢檢驗(yàn)協(xié)調(diào)元雖可以保證總位能從上往下地正確結(jié)果單調(diào)收斂，但往往過于復(fù)雜，使用麻煩。在工程上往往使用形式簡單的非協(xié)調(diào)元，自然，最小位能原理對此不再適用，那么在什么條件下，這類元素才能導(dǎo)致向正確解收斂呢？ Irons 提出了一個稱作“拼片實(shí)驗(yàn)”（patch test） 的檢驗(yàn)方法。實(shí)踐表明，這種檢驗(yàn)方法是有效的，但“拼片實(shí)驗(yàn)”的理論證明尚不清楚。拼片試驗(yàn)內(nèi)容為“假設(shè)

30、由若干元素拼成的一個任意拼片處于等應(yīng)力狀態(tài)，這時，其位移函數(shù)w( x , y ) 一般可用一m階完全多項(xiàng)式函數(shù) P( x , y ) 表示，（入在薄板問題中，m=2） ，而且，在這一拼片的邊界上，也設(shè)置了符合等應(yīng)變狀態(tài)的位移邊界條件。然后，將需要檢驗(yàn)的某種元素按此條件進(jìn)行計(jì)算。如果最后得到的有限元法解答能和P( x , y ) 一致，那么，稱這種元素能夠通過拼片試驗(yàn)，而通過拼片試驗(yàn)的元素將給出收斂的結(jié)果。注：P( x , y ) 至少應(yīng)能代表各種等應(yīng)變狀態(tài)。如： 拼片： 9節(jié)點(diǎn)參數(shù)三角變元 *拼片test實(shí) 際上成為非協(xié)調(diào)元的收斂準(zhǔn)則（1） 注釋： 、兩種拼片均有滿足各種等應(yīng)變狀態(tài)的

31、邊界條件。、用非協(xié)調(diào)元解上述問題。、對比，不一致。、一致，不一致。、說明通過，的網(wǎng)格劃分不通過，但誤差不大，尚可用。、對于更不規(guī)則的網(wǎng)格，其誤差可能更大，大到不宜采用。§ 1、4、3 廣義協(xié)調(diào)元簡介1、 概念： 廣義協(xié)調(diào)元利用對變分原理的改造，使邊界位移在位能條件降低，從而獲得比非協(xié)調(diào)元有更好性質(zhì)的單元，即對于不規(guī)則網(wǎng)格能通過分片檢驗(yàn)；提高計(jì)算精度。2、 原理： 用能量和加權(quán)線量法聯(lián)合改造。、勢能原理： （b）設(shè)單元邊界真實(shí)變量為 ( x , y ),則在單元邊界上每一點(diǎn)若均能滿足，= 0 ( c )則單元協(xié)調(diào)，位移跨越單元連續(xù)，若不能保證（c）式成立，則為非協(xié)調(diào)元，其能量泛函需要做

32、如下修正：) d s ( d )( d ) 式增加了單元不協(xié)調(diào)位移的能量貢獻(xiàn)項(xiàng)。為Lagrange 乘子，其物理意義為單元邊界的函數(shù)。( d )式可通過帶約束的廣義變分原理求解。如果用加權(quán)殘量法，放松單元向位移協(xié)調(diào)條件，即（）式中T為權(quán)函數(shù)，以式（c）為約束條件建立單元位移形函的矩陣，可使單元向位移協(xié)調(diào)條件在某種積分意義下得到滿足。、（e）式實(shí)際上是一個統(tǒng)一的形式，如取T權(quán)函數(shù)為函數(shù)，即得加權(quán)殘方程為：0、全Tn為邊界應(yīng)力，則權(quán)殘方程為：表示常量應(yīng)力，等價于Irons的分片檢驗(yàn)，是非協(xié)調(diào)遠(yuǎn)的收斂準(zhǔn)則。、具體廣義協(xié)調(diào)元的構(gòu)造方法可參見：龍馭球著“新型有限元之引論”§1.5曲線坐標(biāo)系及等

33、參元概念1、 5、0 引在工程有限元分析中，有幾方面的原因需要需要采用高精度單元：、提高計(jì)算精度：、適應(yīng)不規(guī)則的幾何邊界條件；、用較粗的單元數(shù)據(jù)準(zhǔn)備，獲得較高的計(jì)算效益；為獲得高計(jì)算精度采用通常意義下的單元，本身具有某種缺陷：、增加節(jié)點(diǎn)參數(shù)提高位移參數(shù)的階次，使計(jì)算復(fù)雜，并增加求解規(guī)模、只適用于規(guī)則外形的結(jié)構(gòu)，復(fù)雜邊界難以準(zhǔn)確描述。因此，目前多采用等參元在二維問題中得到了極大的成功，以后，在三維問題和板殼計(jì)算中，也取得了很好的效果。舉例說明對受拉帶孔板應(yīng)力集中的計(jì)算結(jié)果（元與等參元對比）單元類型 m單元的相對邊長應(yīng)力集中系數(shù)K 3.02誤差 6810 0.066 0.16 0.008 2.56

34、7 2.806 2.860157.15.36等參元 6 0.0662.8695.06m等分0w等參元的優(yōu)越性：w、精度高（有線性，也有二次等）；、匹配各種復(fù)雜的直線及曲線邊界；、構(gòu)造單元方法統(tǒng)一，推廣其它類型應(yīng)用不困難。4.剛度計(jì)算、公式；K= = =、Guass消元法消去中間節(jié)點(diǎn)，得聚縮剛陣； 、程序系統(tǒng)不再介紹。§2總體剛陣的形成及程序設(shè)計(jì)從有限元的計(jì)算過程可知，形成單剛、集合總剛、置邊界條件、求等效載荷列陣后，才能得到位移，應(yīng)變及應(yīng)力，課程正是按此工作進(jìn)行，最終達(dá)到我們的目的。 §2.1總剛的特征及其組裝時所需要效能的問題1 總剛特性及投放過程總體結(jié)構(gòu)剛陣的形成是各單

35、元剛陣按節(jié)點(diǎn)序排列的集合矩陣（對號入座）。、對稱性：依據(jù) Betti互易原理（講：在結(jié)構(gòu)體的兩點(diǎn)1、2作用單位位移位移，（）先1點(diǎn)作用,再2點(diǎn)作用 ( 2 ) 先2點(diǎn)作用,再1點(diǎn)作用 、帶狀性：（剛陣元素較集中在對角線的位置）一個節(jié)點(diǎn)自由度對應(yīng)總剛的一行（一個節(jié)點(diǎn)不可能和結(jié)構(gòu)的所有節(jié)點(diǎn)相連），這一行總是和相關(guān)的節(jié)點(diǎn)剛陣元素相關(guān)，這也是總剛的組織進(jìn)程反映出的特性。、稀疏性：不相關(guān)的單元（指節(jié)點(diǎn)沒有任何聯(lián)系的單元）沒有共同聯(lián)系的剛性元素，故在對號累加中始終為0。1例： 節(jié)點(diǎn)（ 1 ，4）分屬兩個單元，32 故剛陣中無K系數(shù)，同理 （1，6）、（2，6）、（3，4）、（4，6） 4654 )投放過程

36、：*總剛按節(jié)點(diǎn)排序（不與單元節(jié)點(diǎn)序相關(guān)） 考慮的問題：總的按節(jié)點(diǎn)排序： 2、投放。2、 考慮的問題： 坐標(biāo)轉(zhuǎn)換 ；主從節(jié)點(diǎn)轉(zhuǎn)換；有效方式及其效率，依據(jù)總剛的特點(diǎn)，尋找節(jié)約的有序方法（由于對稱可僅有一半）存貯效率措施使用存貯單元的多少，用總剛規(guī)模替標(biāo)度量?？倓傄?guī)模：存貯總剛單元的總數(shù)及其存貯的形式。*存貯方式 1、三角存放 （一維存） 總剛規(guī)模：KK= 一維存時，要有一個投放算法，使用時，要有一個調(diào)用算法。 n m方式：可一維存，可二維存（不好）切齊即可2維存*存貯方式2、 等帶寬存貯： KK = m x n ；m 為最大帶寬。 ××××× 按節(jié)

37、點(diǎn)編號的最大差計(jì)算；由于帶狀性質(zhì) ×××× m<<n （可一維，也可二維） ×××× ×××× ××××*變帶寬存貯： 可將每行第1個非0元素以外的0元素全部排除在存貯單元之外，需要專門的程序計(jì)算 KK. 只能一維存，還需要一個每行第1非0元素的管理。 *全希疏存放： 僅存貯各單元的非0元素及分解出的非0元素。需要記每一行的首尾特征，行元素的列位置等，則要有復(fù)雜的管理機(jī)構(gòu)。管理機(jī)構(gòu)也消耗單元，效率以遠(yuǎn)小于其它存方式消耗單元總

38、數(shù)為恰當(dāng)。求解效率： 一般不計(jì)編程的效率（程序編的優(yōu)劣不計(jì)），只計(jì)存貯效率和求解效率。求解效率當(dāng)然以機(jī)時計(jì)算 ，效率成本（即存貯求解的消耗成本）。求解效率與存貯效率要綜合評價，但存貯效率一般放在第一位考慮。設(shè)：執(zhí)行的時間為T，存貯量為S，則一個存貯與求解的成本估價是T, S的復(fù)雜函數(shù)，即： C( S, T ) = T * P( S ) (P 是大于1的多項(xiàng)式)。重要的存貯效率，故上述函數(shù)反映了節(jié)省存貯資源的重要性，因?yàn)闇p少存貯要求后，也許不再需要動用外存，從而節(jié)省了數(shù)據(jù)傳輸?shù)馁M(fèi)用。T實(shí)際上也是S的函數(shù)，存入更多的元素，也就需要更多的計(jì)算，因此，該意義上講也要省存貯。3、 影響總剛規(guī)模的因素：

39、除存貯方式的因素之外，另一個影響總剛規(guī)模的因素是節(jié)點(diǎn)編號的順序（按帶狀存貯時，也按全稀疏存貯（編號僅是一個因素），可能影響非0元的分解時的增加策略） 節(jié)點(diǎn)編號對帶寬的影響：變帶寬KK=17總綱形式321 帶寬的主要影響因素是節(jié)點(diǎn)單元間的編號，差越大，則帶寬寬，若差越小，各非0元集中在對角線附近，便會使帶寬減少，這一點(diǎn)構(gòu)成了我們過行節(jié)點(diǎn)編號的基本原則，舉例說明：最大帶寬m = 6KK=36 1 ×2 × ×3 × × ×4 × × × ×5 × × × ×

40、×6 × × × × × ×456m = 4KK=24171 ×2 × ×3 × × ×4 × × × ×5 × × × × ×6 × × ×3216541 ×2 × ×3 × × ×4 × × ×5 × × ×6 ×

41、× ×53115m = 3KK=18246由此可見，取定了簡化模型之后，對節(jié)點(diǎn)如何編號是一個很重要的問題，可直接越小到存貯規(guī)模的大小和機(jī)時的多少，但如何把節(jié)點(diǎn)的編號的順序便的較優(yōu)，是一個很復(fù)雜的問題，這將在以后專門研究這里給出一個編號的原則（工程準(zhǔn)則），即是“優(yōu)先沿結(jié)構(gòu)較寬的方向按坐標(biāo)進(jìn)行編號。 節(jié)點(diǎn)單元編號的自動形成所謂節(jié)點(diǎn)單元編號就是給XE,ME賦值，對于一般結(jié)構(gòu)的分析問題，都有很多的節(jié)點(diǎn)單元，若僅用人力在圖紙上劃分單元和計(jì)算節(jié)點(diǎn)坐標(biāo)，則是很麻煩的，也是易于出錯的，因此，這部分工作盡可能地去讓計(jì)算機(jī)完成。自動編序（又稱自動劃分單元），當(dāng)然也需要把結(jié)構(gòu)的簡化模型取定，然后

42、給某個編序原則，（當(dāng)然應(yīng)盡量靠近編排優(yōu)序的原則，當(dāng)然也可以依靠算法去優(yōu)化），這一原則也需要設(shè)計(jì)者事先在草圖上示意性標(biāo)處，以便檢驗(yàn)這原則的實(shí)施正確性。（也可以用圖顯方式對照草圖檢驗(yàn)），最后即可編制程序。舉一例：yx自動說明自動編序?qū)嵤┓椒ǎ?591317(13)(19)(15)(21)(14)(20)26101418JA(JB=3)(16)(22)37111519(11)(17)(23)(12)(18)(24)48121620JB (分4等)SUBROU XEME( JA , JB , A , B , XE , ME , NP) NP=(JA + 1) * (JB + 1)DIMENDION X

43、E( NP , 1 ) ME( NP , 3 )D 1 N = 1 NPXE( N , 1 ) = B * ( JB - ( N -1 )/( JA + 1 )XE( N ,2 ) = -A * ( JA + 1 ) * (FLOAT ( N - 1 )/( JA + 1) - ( N -1 )/( JA + 1 )D 2 I = 1 JAD 2 J = 1 JBKI = 2 * ( I + ( J - 1 ) * JA ) 產(chǎn)生偶數(shù)單元編號KI = 2 、8、14、20第一行ME ( KI , 1 ) = I + ( J - 1 ) * ( JA + 1 ) + 1ME ( KI , 2 )

44、 = ME( KI , 1 ) - 1ME( KI , 3 ) = ME( KI , 1 ) + JA + 1 ME( KI - 1 , 1 ) = ME ( KI , 1 )ME( KI - 1 , 2 ) = ME ( KI , 3 ) - 1 ME( KI - 1 , 3 ) = ME ( KI , 3 ) RETURN END §2.2 總剛壓縮存貯的實(shí)現(xiàn) 由上述分析可知，變帶寬存貯 KK 最小，顯然最省機(jī)時，但由于帶寬是變的，故只能用一維數(shù)組依次存放總剛陣中各行的第1 個非0元素到對角線的各元素，令該數(shù)組稱為AKM(kk). 要使用AKM(kk)，還必須另外一個整形數(shù)組來記

45、憶原總剛陣中K中每一行（對應(yīng)一個自由度） 的非0 元素到對角線的個數(shù)，該數(shù)組計(jì)為KDKM(NID),NID 為總自由度數(shù)。 在實(shí)用中，我們并不是真正去記憶每行的非0 元素?cái)?shù)，而是記憶每行非0元素的累加數(shù)，也即每行最后一個元素在一維數(shù)組中的位置，這樣給使用帶來了很大的方便。所以在機(jī)器程序中記憶位置的編號，經(jīng)常采用后一種方法。按這種方法，例2的KDKM（I）的元素應(yīng)為1，3，5，9，13，17（KK）。而第3個例題則應(yīng)為1，3，6，9，12，15（KK）。顯然這時KDKM( NID ) = KK (最后一個數(shù))。1 計(jì)算原理： P， 點(diǎn)。解釋：行快是以每一個節(jié)點(diǎn)所對應(yīng)的自由度數(shù)為行數(shù)，列號從1到最

46、后一列，2維是 Row = 2 ; 3 維是 Row = 3 .A、 只有與節(jié)點(diǎn)之構(gòu)成的相關(guān)點(diǎn)（與節(jié)點(diǎn)相關(guān)的那些單元的編號），才在總剛 陣中的第i 行塊上提供非0子塊。B、 第i行塊上非0子塊的分布情況與節(jié)點(diǎn) i的相關(guān)節(jié)點(diǎn)編號的差值有關(guān)。帶寬對于下 陣，則是相關(guān)節(jié)點(diǎn)最小編號列對角線子塊。對上陣存貯，則是相關(guān)節(jié)點(diǎn)中最大編號列對角線子塊（從右到左）。2、框圖：書P 。2. 出口剛陣的形成計(jì)算由式再由Gauss數(shù)值積分，可知幾個分塊剛陣的計(jì)算如下： Kd已在前節(jié)介紹過了。于是取代計(jì)算為：計(jì)算出后，即可按式求解出口剛陣上述是理論求解過程，要求的逆及矩陣乘，矩陣求解是費(fèi)時的，因此，一般真實(shí)計(jì)算時，并不沿循理論理論求解過程而是采用高斯消去方法，下面講述此過程。由式可知，為單元內(nèi)節(jié)點(diǎn)的等效載荷列陣（實(shí)際并不指出是那個具體節(jié)點(diǎn)，相當(dāng)于不增加節(jié)點(diǎn)。），因形成時僅對單元邊界節(jié)點(diǎn)進(jìn)行，故，所以式可寫為： 利用矩陣的初等變換，消去陣，使在的