非線性有限元及彈塑性力學(xué)chap1-典尚設(shè)計(jì)

1、第一章 非線性代數(shù)方程組的數(shù)值解法1.1 直接迭代法1.2 牛頓法和修正牛頓法1.3 擬牛頓法1.4 增量方法1.5 增量弧長(zhǎng)法2000.3哈爾濱工業(yè)大學(xué) 王煥定教授制作2 非線性問題可分為三類：材料非線性 不管那類非線性問題，最終都?xì)w結(jié)為一組非線性方程(a)=0，a為待求的未知量。 對(duì)許多問題，用某些方法可將(a)=0改造成(a) =P(a)-R=K(a) a -R=0的形式。 對(duì)非線性問題的方程(a)=0，一般只能用數(shù)值方法求近似解答。 、幾何非線性 和邊界非線性。 我們只討論前兩類問題。 其實(shí)質(zhì)是，用一系列線性方程組的解去逼近所討論非線性方程組的解。本章將簡(jiǎn)單介紹有限元分析中常見的各種求

2、解非線性方程組的數(shù)值方法。2000.3哈爾濱工業(yè)大學(xué) 王煥定教授制作31.1 直接迭代法 當(dāng)用某些方法將(a)=0改造成迭代格式(a) =P(a)-R=K(a) a -R=0后a1= K(a0)-1R如果問題是收斂的， a1將比a0有所改善。an+1= K(an)-1R an=an+1- an當(dāng)設(shè)范數(shù)為或設(shè)范數(shù)為收斂條件則為 ，設(shè)一初始未知量a0 ，則由它可得 如此反復(fù)迭代可得2000.3哈爾濱工業(yè)大學(xué) 王煥定教授制作4 如果考慮到每步迭代(an) =P(an)-R=K(an) an -R0將(an)視為不平衡力（或失衡力）并作為衡量收斂的標(biāo)準(zhǔn) 應(yīng)指出的是，對(duì)單變量情況，如講義圖示，直接迭代實(shí)

3、質(zhì)是“割線”法1.1 直接迭代法返首頁(yè) ，則收斂條件也可改為 ，一定條件下這種迭代過程是收斂的 ，但對(duì)多自由度情況，由于未知量通過矩陣K(an)的元素互相耦合，在迭代過程中往往出現(xiàn)不穩(wěn)定現(xiàn)象。2000.3哈爾濱工業(yè)大學(xué) 王煥定教授制作51.2 牛頓法和修正牛頓法 如果將非線性方程(a) =0在an 附近展開，則又如果(a)n的逆存在，則an 近似等于 記 KT(an)=(a)n，Pn =(an)an-(a)n-1(an)則 an-KT(an)-1 Pn ， an+1=an+an切線矩陣不平衡力 如此逐步計(jì)算，即可得到非線性方程的解答，這就是牛頓-拉夫森法。(a) =(an)+ (a)n an+

4、。 =0或用求和約定可寫為2000.3哈爾濱工業(yè)大學(xué) 王煥定教授制作61.2 牛頓法和修正牛頓法 牛頓法要每步都計(jì)算切線矩陣KT（也稱剛度）并解線性方程組，雖精度高，但工作量也大。其中n的作用是改變切線矩陣KT的主對(duì)角元素，使奇異性或病態(tài)得到改善。更多的改進(jìn)方法可參看沈聚敏鋼筋混凝土有限元與板殼極限分析等。 此外，在某些非線性問題(如理想塑性和軟化塑性問題)中用牛頓法，迭代過程中切線矩陣可能是奇異的或病態(tài)的 ，為了克服這一現(xiàn)象，可有多種處理方法，其一是按下式來(lái)求2000.3哈爾濱工業(yè)大學(xué) 王煥定教授制作71.2 牛頓法和修正牛頓法 如果在計(jì)算的每一步內(nèi)，矩陣KT都用初始近似解KT0計(jì)算，在這種

5、情況下，僅第一步迭代需要完全求解一個(gè)線性方程組，如果將KT0三角分解并存儲(chǔ)起來(lái)，而以后各步迭代中采用公式則只需對(duì)上式右端項(xiàng)中的 進(jìn)行回代就行了。這種方法稱為修正的牛頓法。 為了提高修正牛頓法的收斂速度可采用某些過量修正技術(shù)。講義上作了簡(jiǎn)要介紹，請(qǐng)大家自己看。返首頁(yè)2000.3哈爾濱工業(yè)大學(xué) 王煥定教授制作81.3 擬牛頓法擬牛頓法的主要思想是： 首先設(shè)(KT)n+1可寫成如下修正形式 接著設(shè)(KT)n+1必須滿足如下所謂擬牛頓方程(KT)n+1=(KT)n+(KT)n式中(KT)n稱為修正矩陣。由此可建立擬牛頓法迭代格式(略去了下標(biāo)T)2000.3哈爾濱工業(yè)大學(xué) 王煥定教授制作9 要用擬牛頓法

6、，還需給出修正矩陣的計(jì)算。推導(dǎo)修正矩陣算式的思路是：(un)和(vn)是秩1（或秩2，講義為秩2）的列向量，將修正矩陣代入擬牛頓方程可得 設(shè) (KT)n=(un)(vn)T如果取(vn)=(a)n，則當(dāng)(a)n(0)時(shí)(KT)n +(un)(vn)T(a)n=()n(un)=(vn)T(a)n-1()n-(KT)n(a)n假設(shè)(vn)T(a)n0，則有(KT)n =(a)n)T(a)n-1()n-(KT)n(a)n (a n)T當(dāng)(a)n=(0)時(shí)，迭代已收斂，(KT)n =(0)。2000.3哈爾濱工業(yè)大學(xué) 王煥定教授制作10 講義上的內(nèi)容比這里說明的多，但基本思路是一樣的。關(guān)于秩2的算法，

7、請(qǐng)大家自己看。 1.設(shè)(a)0求(KT)0 ； 對(duì)秩1算法來(lái)說，實(shí)際使用的步驟為： 2.求 (a)0=-(KT)0-1()0；a1=a0+a0 5.計(jì)算(KT)1=(KT)0+(KT)0 ;(KT)0= (KT)1. 3. 計(jì)算()0； a0=a1 。 6.重復(fù)第2步，直到達(dá)到精度要求為止。 4. 計(jì)算 (KT)0 = ()0-(KT)0(a)0 (a0)T/(a0)T(a)0； 講義上的塞爾曼公式可用逆矩陣定義驗(yàn)證。2000.3哈爾濱工業(yè)大學(xué) 王煥定教授制作11返首頁(yè)講義上給出了三種方法的對(duì)比，指出了選用什么算法應(yīng)考慮的因素。2000.3哈爾濱工業(yè)大學(xué) 王煥定教授制作12 求解非線性方程組的

8、另一類方法是增量方法。 當(dāng)問題的性質(zhì)與加載的歷史有關(guān)時(shí)，如彈塑性問題，則必須采用增量方法。1.4 增量方法使用這種方法需要知道“荷載”項(xiàng)(R)為零時(shí)問題的解(a)0。 在實(shí)際問題中，(R)經(jīng)常代表真實(shí)荷載，(a)0 代表結(jié)構(gòu)位移。在問題的初始狀態(tài)，它們均為零。 這種從問題的初值開始，隨著荷載列陣(R)按增量形式逐漸增大，研究(a)i的變化規(guī)律的方法，稱為增量方法。2000.3哈爾濱工業(yè)大學(xué) 王煥定教授制作13 設(shè)“荷載”(R)在任一增量步的值為(R)， 為荷載增量因子，(R)為標(biāo)準(zhǔn)荷載列陣,則非線性方程 (a)= (0)可寫為 引入切線矩陣且略去高階小量后可改寫為(a,)=(P(a)-(R)=

9、(0)若+時(shí)的解答為(a)+(a)，象牛頓法一樣，將( (a)+(a),+) 按Taylor級(jí)數(shù)展開，則可得2000.3哈爾濱工業(yè)大學(xué) 王煥定教授制作14 設(shè)荷載增量因子分別取如下值(a)m+1=(a)m+(a)m則荷載(R)可分成M級(jí)，第m級(jí)荷載為m(R)，其增量為(m+1-m)(R)=m(R)。 由此可得(a)m=KT(a)m,m)-1m(R)但是，這樣做的每一步都將產(chǎn)生誤差，結(jié)果使解答漂移。講義上簡(jiǎn)單介紹了四種解決漂移的方法，下面僅對(duì)混合法作簡(jiǎn)單說明，其他方法請(qǐng)大家自行閱讀。2000.3哈爾濱工業(yè)大學(xué) 王煥定教授制作15 所謂混合法是指，在增量法每一增量步進(jìn)行自修正的迭代計(jì)算。其m增量步

10、n次迭代的計(jì)算公式為 在實(shí)際計(jì)算中，對(duì)于 mM-1的各增量步的計(jì)算，可以只進(jìn)行少許幾次(例如3次)迭代，而對(duì)于m=M-1，即最后的一個(gè)荷載增量，需耍使用較多次迭代，以使近似解更接近于真解。返首頁(yè)自修正不平衡力 用混合法求解時(shí)，所選取的荷載增量的步長(zhǎng)可以比普通增量算法的步長(zhǎng)大一些。2000.3哈爾濱工業(yè)大學(xué) 王煥定教授制作16 用迭代法或增量法進(jìn)行極限分析時(shí)，在極值點(diǎn)附近往往可能不收斂。這時(shí)可用增量弧長(zhǎng)法來(lái)解決。 為便于理解，以桿單向拉伸為例加以說明。1.5 增量弧長(zhǎng)法 增量弧長(zhǎng)法的基本思想是：將作為獨(dú)立變量，在每個(gè)增量步進(jìn)行自修正法平衡迭代，在迭代過程中自動(dòng)控制荷載因子的取值。也即前步結(jié)果本步

11、n次增量2000.3哈爾濱工業(yè)大學(xué) 王煥定教授制作17 如圖所示，矢徑可表達(dá)為uuu因?yàn)?由于弧長(zhǎng)法引入了如下約束方程因此由此可得2000.3哈爾濱工業(yè)大學(xué) 王煥定教授制作18 由矢量代數(shù)和約束方程可得因此若記也即則2000.3哈爾濱工業(yè)大學(xué) 王煥定教授制作19將其代入約束方程，可得式中系數(shù)為上述式子是從簡(jiǎn)單情況推出的，如果除 外均理解為矩陣，即為一般情況的弧長(zhǎng)法方程。 一元二次方程有兩個(gè)根，應(yīng)取 和 間成銳角的根。據(jù)此可建立判別條件，具體推導(dǎo)這里從略了。2000.3哈爾濱工業(yè)大學(xué) 王煥定教授制作20這樣做迭代的軌跡很接近圓弧，而計(jì)算工作量可減少很多。 從上述公式可見，求 的工作量是很大的，為此，可令 和 相互垂直，也即2000.3哈爾濱工業(yè)大學(xué) 王煥定教授制作21 由相互垂直的條件可得 綜上所述，弧長(zhǎng)法求解步驟為：1)選定荷載參考值 ，和本步荷載因子 ，解得 ，由 求弧長(zhǎng)。2)修改切線剛度矩陣并三角化。檢查對(duì)角元，正定則加載，負(fù)定則加負(fù)荷載。若矩陣對(duì)應(yīng)行列式為零，達(dá)到極限荷載。4)由 和 求 和