基于位移控制新法的結(jié)構(gòu)非線性有限元數(shù)值模擬

陳夢(mèng)成，李嘉鈺

（華東交通大學(xué)土木建筑學(xué)院,江西 南昌330013）

結(jié)構(gòu)全過程分析是非線性分析理論的基礎(chǔ)，它常常應(yīng)用于大跨度空間結(jié)構(gòu)的設(shè)計(jì)之中。 在結(jié)構(gòu)非線性全過程分析中，剛度矩陣的組裝、計(jì)算方法及迭代方式的選取，都將對(duì)結(jié)果產(chǎn)生重大的影響。 近年來，人們一直在對(duì)結(jié)構(gòu)非線性全過程分析中的迭代控制算法進(jìn)行研究，提出了壓縮平衡迭代法、人工彈簧法、位移控制法、弧長(zhǎng)法等多種方法。 王志軍等[1]在常規(guī)的位移控制法和弧長(zhǎng)法中引入了簡(jiǎn)單適用的考慮路徑搜尋和單因子加速的方法，能明顯改善收斂情況；鄭宏等[2]利用Sherman-Morrison 定理，對(duì)經(jīng)典的荷載增量控制法程序做了少量的修補(bǔ)，構(gòu)建了一個(gè)有良好數(shù)值特性且僅需很少運(yùn)算量的位移控制法程序；杜柏松等[3]引入廣義位移作為約束方程來確定荷載增量因子，從而控制結(jié)構(gòu)荷載加載的全過程，形成了廣義的位移控制法；許茂等[4]針對(duì)負(fù)剛度結(jié)構(gòu)，采用了附加正剛度彈簧的方法進(jìn)行了分析，解決了負(fù)剛度切線矩陣在有限元分析中造成的迭代求解不收斂的問題；李鋼等[5]證明了擬力法在計(jì)算結(jié)構(gòu)負(fù)剛度時(shí)的優(yōu)勢(shì)，它是一種高效的非線性分析方法；曲激婷等[6]基于向量式有限元的基本思想，提出了力-位移控制的復(fù)合模擬方法，此方法可以很好地求解結(jié)構(gòu)的大變形問題；AL－Aukaily 等[7]采用了直接微分法和向前差分法，同時(shí)考慮結(jié)構(gòu)的幾何非線性和材料非線性，對(duì)求解結(jié)構(gòu)響應(yīng)所用位移控制法的相關(guān)參數(shù)進(jìn)行了敏感度分析。

利用杜修力等提出的位移控制新方法[8]，編寫相關(guān)計(jì)算程序，并結(jié)合實(shí)際算例，對(duì)結(jié)構(gòu)的非線性全過程進(jìn)行有限元分析，驗(yàn)證位移控制新方法的有效性、正確性和可行性。

1 U.L.列式下增量有限元方程

根據(jù)連續(xù)性介質(zhì)力學(xué)和虛位移原理[9]，并以時(shí)刻結(jié)構(gòu)的狀態(tài)作為基準(zhǔn)參考構(gòu)型，應(yīng)變、應(yīng)力以及位移的增量關(guān)系為[10]

式中：KL為整體坐標(biāo)下的彈性剛度矩陣；KG為整體坐標(biāo)下的幾何剛度矩陣，兩者之和為單元切線剛度矩陣；fe為單元抵抗力向量。

2 位移控制新方法

2.1 傳統(tǒng)位移控制法

傳統(tǒng)的位移控制法（displacement control method，DCM）是求解材料、構(gòu)件及結(jié)構(gòu)包含極限強(qiáng)度后軟化下降段的荷載-位移關(guān)系非線性全過程曲線的常用方法之一，它是通過引入給定力作用節(jié)點(diǎn)（位移控制點(diǎn)）的位移增量來反求荷載增量和其它未知節(jié)點(diǎn)的位移增量。 其一般迭代方程為[11]

式中：[KT]為結(jié)構(gòu)的切線剛度矩陣；Δλ 為荷載因子增量；{P}為參考荷載向量；[R]為不平衡力向量；j 為增量步數(shù)；i 為迭代次數(shù)。 接下來根據(jù)需要控制的位移向量對(duì)結(jié)構(gòu)的迭代平衡方程式（22）進(jìn)行重新排序，即對(duì)結(jié)構(gòu)的切線剛度矩陣、荷載向量以及不平衡力向量進(jìn)行如下分塊

至此，可以計(jì)算出每一步迭代中的荷載因子增量Δλ，最終可計(jì)算其余點(diǎn)位移分量。

2.2 位移控制新方法

傳統(tǒng)位移控制法的求解過程需將有限元方程的剛度矩陣進(jìn)行分塊并重新排列，從而得到新的方程和系數(shù)矩陣，但新的系數(shù)矩陣將不再具有原剛度矩陣的對(duì)稱性及帶狀性，會(huì)使計(jì)算量顯著增加，且計(jì)算過程較為復(fù)雜，最終帶來諸多計(jì)算上的不便。

一般地，有限元增量方程可寫為

3 算例分析

3.1 幾何非線性

選取William 肘式框架[12]進(jìn)行計(jì)算分析，William 肘式框架是眾多學(xué)者檢驗(yàn)其所建立理論有效性及正確性的典型算例，其原因在于該結(jié)構(gòu)具有顯著的幾何非線性特征，同時(shí)又具備可以用來與理論值對(duì)比的試驗(yàn)結(jié)果，如圖1 所示。 該結(jié)構(gòu)兩端固結(jié)，水平方向跨度為657.6 mm， 豎直方向高為9.8 mm，桿件為細(xì)長(zhǎng)矩形等截面直桿，截面寬19.1 mm，高6.17 mm，彈性模量為71 018.5 MPa，框架頂端作用一豎向集中力荷載。

在分析時(shí)， 將該結(jié)構(gòu)離散成20 個(gè)Euler-Bernoulli[13]梁?jiǎn)卧M(jìn)行模擬，取控制點(diǎn)的位移增量步步長(zhǎng)u=0.1 mm，總位移增量步n=180 即分析結(jié)構(gòu)位移控制點(diǎn)0 至18 mm 的荷載-位移關(guān)系曲線，允許位移收斂誤差δ=10-4，在U.L.列式下，使用MATLAB 進(jìn)行編程計(jì)算。

圖2 給出了運(yùn)用傳統(tǒng)位移控制法、 位移控制新方法得到的數(shù)值分析結(jié)果與William 試驗(yàn)結(jié)果之間的對(duì)比；圖3 給出了運(yùn)用傳統(tǒng)位移控制法與位移控制新方法得到的結(jié)構(gòu)最終變形曲線?？梢姡灰瓶刂菩路椒芸缭胶奢d-位移關(guān)系中的極值點(diǎn)以及結(jié)構(gòu)極限強(qiáng)度后的軟化下降段模擬，且能準(zhǔn)確地計(jì)算出結(jié)構(gòu)的荷載-位移全過程曲線。兩種方法所得到的結(jié)構(gòu)最終變形曲線高度吻合，證明了位移控制新方法的可行性和正確性。

圖2 William 肘式框架荷載-位移曲線Fig.2 The load-displacement curve for William frame

圖3 William 肘式框架變形曲線Fig.3 The deformation curve for William frame

3.2 材料非線性

一混凝土棱柱體試件[14]尺寸為150 mm×150 mm×300 mm，混凝土等級(jí)為C30，考慮材料非線性，對(duì)該構(gòu)件進(jìn)行非線性有限元分析。

計(jì)算分析時(shí)，采用三維纖維梁?jiǎn)卧?，沿長(zhǎng)度方向離散成3 個(gè)單元， 截面劃分成15×15 的纖維網(wǎng)格，取控制點(diǎn)的位移增量步步長(zhǎng)u=0.1 mm，總位移增量步n=24，即分析構(gòu)件應(yīng)變?yōu)?～0.008 時(shí)對(duì)應(yīng)的應(yīng)力大小，允許位移收斂誤差δ=10-4?；炷恋谋緲?gòu)模型根據(jù)我國(guó)混凝土設(shè)計(jì)規(guī)范[15]，采用單軸受壓應(yīng)力-應(yīng)變曲線。

圖4 混凝土柱應(yīng)力-應(yīng)變曲線Fig.4 Stress-strain curve of concrete column

將運(yùn)用傳統(tǒng)位移控制法和位移控制新方法得到的數(shù)值分析結(jié)果與文獻(xiàn)[14]結(jié)果進(jìn)行對(duì)比，如圖4 所示?？梢姡瑑煞N方法計(jì)算出的曲線與文獻(xiàn)中曲線的吻合程度良好，且兩種方法的結(jié)果高度吻合，證明了位移控制新方法的準(zhǔn)確性和可行性；除此之外，本文在利用MATLAB 進(jìn)行編程計(jì)算時(shí)，對(duì)兩種方法的計(jì)算完成時(shí)間進(jìn)行了測(cè)試，結(jié)果發(fā)現(xiàn)，傳統(tǒng)位移控制法所需時(shí)間約為333 s，而位移控制新方法所需時(shí)間約為202 s，新方法近節(jié)約了40%的時(shí)間，優(yōu)化了計(jì)算過程。

4 結(jié)論

通過以上理論分析與數(shù)值模擬計(jì)算，可以得到以下結(jié)論：

1） 位移控制法可以模擬出結(jié)構(gòu)的荷載-位移全過程非線性曲線，尤其是可以模擬出曲線含有極值點(diǎn)或軟化下降段的情況；

2） 位移控制新方法不會(huì)破壞有限元?jiǎng)偠染仃嚨膶?duì)稱性和帶狀性，且能得到穩(wěn)定、精確的數(shù)值解；

3） 位移控制新方法較傳統(tǒng)位移控制法更為簡(jiǎn)便實(shí)用，求解所需的計(jì)算步驟和計(jì)算時(shí)間都可以得到一定程度上的優(yōu)化，節(jié)約了計(jì)算成本。