孫寶印，申 偉，孫天舒，歐進(jìn)萍,4

(1. 河海大學(xué)土木與交通學(xué)院，江蘇，南京 210098；2. 廣西大學(xué)土木建筑工程學(xué)院，廣西，南寧 530004；3. 大連理工大學(xué)土木工程學(xué)院，遼寧，大連 116024；4. 哈爾濱工業(yè)大學(xué)土木工程學(xué)院，黑龍江，哈爾濱 150090)

計(jì)算效率和模擬精度是高層建筑等重大工程地震彈塑性全過(guò)程分析的關(guān)鍵技術(shù)瓶頸。一方面，此類(lèi)結(jié)構(gòu)的規(guī)模龐大，為了精確刻畫(huà)其非線性行為，有必要針對(duì)整體結(jié)構(gòu)建立精細(xì)化的有限元模型，自由度動(dòng)輒數(shù)十上百萬(wàn)，結(jié)構(gòu)剛度矩陣的集成和分解耗時(shí)嚴(yán)重。更為關(guān)鍵的是，結(jié)構(gòu)一旦進(jìn)入非線性狀態(tài)，傳統(tǒng)的變剛度法需要反復(fù)迭代更新結(jié)構(gòu)剛度，計(jì)算效率更是成為棘手難題。另一方面，高層結(jié)構(gòu)地震倒塌破壞主要是由局部關(guān)鍵構(gòu)件嚴(yán)重?fù)p傷導(dǎo)致，而這部分關(guān)鍵構(gòu)件的非線性變形特征具有顯著的多樣性和復(fù)雜性，對(duì)模型精細(xì)程度提出了較高的要求。實(shí)際工程中為了滿足計(jì)算效率要求，無(wú)法針對(duì)整體結(jié)構(gòu)建立精細(xì)模型，進(jìn)而難以保證局部關(guān)鍵構(gòu)件的模擬精度。由此可見(jiàn)，計(jì)算效率和模擬精度之間相互矛盾、制約?，F(xiàn)有的計(jì)算方法較難有效兼顧這類(lèi)大規(guī)模結(jié)構(gòu)非線性分析的效率和精度，因此，發(fā)展更加高效、精確的非線性分析理論和計(jì)算方法十分必要。

高層結(jié)構(gòu)在地震作用下呈現(xiàn)局部非線性特征，如何充分利用這種局部化特征來(lái)提高大規(guī)模結(jié)構(gòu)非線性分析的計(jì)算效率，是近幾十年來(lái)國(guó)內(nèi)外學(xué)者關(guān)注的重點(diǎn)。Clough 等[1]利用子結(jié)構(gòu)方法對(duì)局部非線性的大規(guī)模結(jié)構(gòu)系統(tǒng)進(jìn)行動(dòng)力分析。將整體結(jié)構(gòu)劃分為線彈性和塑性區(qū)域：一方面，通過(guò)靜力凝聚消去線彈性區(qū)域內(nèi)部自由度，降低求解矩陣規(guī)模來(lái)提高計(jì)算效率；另一方面，提高塑性區(qū)模型精細(xì)程度來(lái)有效兼顧模擬精度。Wong 等[2]提出擬力法，定義能夠描述結(jié)構(gòu)構(gòu)件非線性行為的塑性自由度，建立相應(yīng)的塑性機(jī)制，來(lái)模擬結(jié)構(gòu)的非線性行為。并采用初始剛度對(duì)狀態(tài)方程進(jìn)行迭代求解，利用顯式求解算法，可實(shí)現(xiàn)高效、穩(wěn)定的結(jié)構(gòu)非線性分析。隨后，Li 等[3-5]在擬力法的基礎(chǔ)上作了進(jìn)一步改進(jìn)，將材料應(yīng)變分解為線彈性應(yīng)變和非線性應(yīng)變兩部分，并在單元內(nèi)構(gòu)造非線性應(yīng)變的分布場(chǎng)模型，提出了單元非線性變形機(jī)制的模擬方法。再結(jié)合虛功原理，建立具有隔離非線性特征的結(jié)構(gòu)控制方程，在該方程迭代分析時(shí)僅需對(duì)一個(gè)小規(guī)模的Schur 補(bǔ)矩陣進(jìn)行迭代更新，從而能夠顯著提高結(jié)構(gòu)非線性分析的計(jì)算效率。此外，陸新征等[6]、Sun 等[7-8]以及陳宇等[9]從多尺度建模角度，對(duì)整體結(jié)構(gòu)進(jìn)行差異化建模分析，在保證模擬精度的前提下，可大幅降低結(jié)構(gòu)模型自由度。

傳統(tǒng)的子結(jié)構(gòu)方法、擬力法以及多尺度方法均需要預(yù)判結(jié)構(gòu)的塑性區(qū)域，并且模型網(wǎng)格一旦確認(rèn)，分析過(guò)程中較難更改。為了規(guī)避這些缺陷，孫寶印等[10-13]提出數(shù)值子結(jié)構(gòu)方法(如圖1所示)，將原始高層結(jié)構(gòu)的非線性分析轉(zhuǎn)化為主結(jié)構(gòu)的等效線彈性分析和屈服構(gòu)件的隔離子結(jié)構(gòu)非線性分析。在整個(gè)地震彈塑性時(shí)程分析過(guò)程中，由線彈性主結(jié)構(gòu)模型計(jì)算得到整體結(jié)構(gòu)的動(dòng)力響應(yīng)，根據(jù)結(jié)構(gòu)全局位移實(shí)時(shí)判斷構(gòu)件的屈服狀態(tài)。將構(gòu)件屈服隔離并建立精細(xì)化的子結(jié)構(gòu)模型，由子結(jié)構(gòu)非線性分析計(jì)算得到屈服構(gòu)件的真實(shí)內(nèi)力，并將其傳回主結(jié)構(gòu)。主結(jié)構(gòu)將轉(zhuǎn)化得到的屈服構(gòu)件等效外荷載作為外力項(xiàng)施加在屈服構(gòu)件端部節(jié)點(diǎn)上，最后，再由主結(jié)構(gòu)進(jìn)行等效線彈性分析得到整體結(jié)構(gòu)的真實(shí)響應(yīng)。經(jīng)研究表明：該方法可以有效兼顧大規(guī)模結(jié)構(gòu)非線性分析的計(jì)算效率和模擬精度。

由于上述數(shù)值子結(jié)構(gòu)方法采用初始彈性剛度進(jìn)行非線性迭代分析，仍存在一定不足。于是，本文擬對(duì)該方法作進(jìn)一步改進(jìn)，提出一種降階的牛頓迭代數(shù)值子結(jié)構(gòu)方法。首先，簡(jiǎn)介數(shù)值子結(jié)構(gòu)方法基本原理及其不足；其次，提出降階的牛頓迭代算法，并介紹屈服構(gòu)件的隔離子結(jié)構(gòu)非線性分析；最后，通過(guò)一平面15 層3 跨鋼結(jié)構(gòu)的地震彈塑性時(shí)程分析來(lái)驗(yàn)證所提方法的精度和效率。

1 數(shù)值子結(jié)構(gòu)方法基本原理及不足

1.1 數(shù)值子結(jié)構(gòu)方法

根據(jù)數(shù)值子結(jié)構(gòu)方法基本原理[10-13]，結(jié)構(gòu)運(yùn)動(dòng)方程可以表示為：

式中：U、和分別為結(jié)構(gòu)的位移、速度和加速度；M、C0、K0、F和Fc分別為結(jié)構(gòu)的質(zhì)量矩陣、初始阻尼矩陣、初始剛度矩陣、外力和非線性修正力，分別表示為：

式中：Ne和Ns分別為總單元數(shù)和屈服單元數(shù)；B(e)為用于檢索結(jié)構(gòu)中單元所在位置的布爾矩陣；為屈服單元的非線性修正力，可表示為：

式中：ke和ue分別為屈服單元的初始彈性剛度和節(jié)點(diǎn)位移；re為由隔離子結(jié)構(gòu)計(jì)算得到的屈服單元內(nèi)力。

式(1)等式左邊是一個(gè)線彈性系統(tǒng)，結(jié)構(gòu)非線性是由等式右邊的非線性修正力Fc來(lái)考慮。對(duì)于局部非線性結(jié)構(gòu)系統(tǒng)，屈服單元數(shù)Ns＜＜Ne遠(yuǎn)小于總單元數(shù)Ne，則非線性修正力Fc的計(jì)算量較小。

1.2 數(shù)值子結(jié)構(gòu)方法的不足

以文獻(xiàn)[14]中15 層3 跨的平面鋼框架結(jié)構(gòu)為研究對(duì)象。采用不動(dòng)點(diǎn)迭代數(shù)值子結(jié)構(gòu)方法[10-13]對(duì)其進(jìn)行地震彈塑性時(shí)程分析時(shí)，主結(jié)構(gòu)中每根梁柱構(gòu)件由30 個(gè)基于歐拉-伯努利梁理論的彈性單元模擬，分析過(guò)程中根據(jù)截面最外側(cè)應(yīng)變是否超過(guò)材料屈服應(yīng)變來(lái)判斷構(gòu)件的彈塑性狀態(tài)[15]，屈服單元由隔離子結(jié)構(gòu)進(jìn)行非線性分析。在隔離子結(jié)構(gòu)中，屈服單元采用基于歐拉-伯努利梁理論的分布式塑性鉸纖維單元模擬，鋼材采用OpenSees中Steel02 本構(gòu)模型模擬，材料彈性模量E=2.0×106MPa，硬化比b=0.1。在不同工況和收斂誤差下，不動(dòng)點(diǎn)迭代數(shù)值子結(jié)構(gòu)方法的迭代次數(shù)時(shí)程曲線如圖2 所示，迭代信息如表1 所示?？梢钥闯?，當(dāng)收斂誤差為1.0×10-5mm 時(shí)，8 度和9 度罕遇地震作用下，數(shù)值子結(jié)構(gòu)方法需要最多24 步和26 步收斂，平均迭代次數(shù)分別為4.6 和6.9；而當(dāng)收斂誤差減小至1.0×10-8mm 時(shí)，最多需要41 步和44 步，平均迭代次數(shù)分別為7.1 和11.1。當(dāng)荷載水平增大或收斂誤差減小，不動(dòng)點(diǎn)迭代數(shù)值子結(jié)構(gòu)方法所需的迭代收斂次數(shù)均有所增加，尤其是當(dāng)收斂條件變得苛刻，收斂速率顯著降低。于是，下文將對(duì)該方法作進(jìn)一步改進(jìn)，以提高其收斂速率。

圖2 不動(dòng)點(diǎn)迭代數(shù)值子結(jié)構(gòu)方法迭代次數(shù)時(shí)程曲線Fig. 2 Iteration number time-history curves using the numerical substructure method based on fixed-point iteration

表1 不動(dòng)點(diǎn)迭代數(shù)值子結(jié)構(gòu)方法迭代信息統(tǒng)計(jì)表Table 1 Table of iteration information using the numerical substructure method based on fixed-point iteration

2 降階的牛頓迭代算法

經(jīng)過(guò)時(shí)間離散后，第n時(shí)步結(jié)構(gòu)的等效運(yùn)動(dòng)方程滿足：

不失一般性，以Newmark-β 逐步積分法為例，第n時(shí)步速度和位移分別表示為：

式中， γ 和 β為控制參數(shù)。本文采用平均加速度法，則 γ 和 β取值分別為0.5 和0.25。

由式(5)可以得到第n時(shí)步速度和加速度：

式 中：c1=γ/(βΔt)；c2=1/β(Δt)2；a1=1-γ/β；a2=(1-γ/2β)Δt；a3=-1/βΔt；a4=1-1/2β。

將式(6)代入式(4)，得到結(jié)構(gòu)的等效平衡方程：

式中，K和Pn分別為結(jié)構(gòu)等效初始剛度和外力：

由式(7)可以得到結(jié)構(gòu)的位移場(chǎng)：

定義與屈服單元有關(guān)的自由度為塑性自由度，則塑性自由度位移us,n可以表示為：

式中，Bs,n為檢索塑性自由度位置的布爾矩陣。

將式(10)代入式(11)，塑性自由度位移殘差φ(us,n)可以表示為：

采用牛頓迭代算法，已知當(dāng)前塑性自由度位移us,n和殘差φ(us,n)，得到更新后的塑性自由度位移：

式中，Hs,n為與塑性自由度有關(guān)的迭代矩陣：

3 隔離子結(jié)構(gòu)非線性分析

考慮到高層結(jié)構(gòu)在地震作用下發(fā)生倒塌破壞，局部關(guān)鍵構(gòu)件非線性程度較高，而大部分屈服構(gòu)件的非線性程度較低(稱之為一般屈服構(gòu)件)。為了滿足模擬精度需求，一般屈服構(gòu)件和少數(shù)關(guān)鍵構(gòu)件所采用的隔離子結(jié)構(gòu)模型的精細(xì)程度不同。為了簡(jiǎn)化起見(jiàn)，本節(jié)以梁柱構(gòu)件為例說(shuō)明。

3.1 一般屈服梁柱構(gòu)件

針對(duì)一般屈服梁柱構(gòu)件，纖維模型可以較為準(zhǔn)確模擬其非線性行為。由于主、子結(jié)構(gòu)中梁柱構(gòu)件均為桿系模型，它們的邊界自由度一致，因此，由隔離子結(jié)構(gòu)可直接計(jì)算得到屈服構(gòu)件內(nèi)力和切線剛度，這里不作詳細(xì)贅述。

3.2 關(guān)鍵梁柱構(gòu)件

以如圖3(a)所示的平面框架結(jié)構(gòu)為例，將主結(jié)構(gòu)中關(guān)鍵柱構(gòu)件隔離并建立精細(xì)化的子結(jié)構(gòu)模型(例如：實(shí)體單元)，如圖3(b)所示。為了消除剛體位移的影響，將主結(jié)構(gòu)中屈服構(gòu)件對(duì)應(yīng)的彈性單元節(jié)點(diǎn)位移ue轉(zhuǎn)換得到局部坐標(biāo)位移ve：

圖3 關(guān)鍵構(gòu)件隔離子結(jié)構(gòu)精細(xì)化模型邊界處理Fig. 3 Boundary treatment of isolated refined substructures for key components

式中：上標(biāo)i為精細(xì)化子結(jié)構(gòu)第i個(gè)邊界；轉(zhuǎn)換矩陣T可以表示為：

式中， α為桿系單元軸向與全局坐標(biāo)系X軸之間的夾角。

如圖3(c)所示，在局部坐標(biāo)系下，根據(jù)平截面假定，精細(xì)化子結(jié)構(gòu)中第i個(gè)邊界中第j個(gè)節(jié)點(diǎn)的位移滿足以下轉(zhuǎn)換關(guān)系：

將式(16)代入式(18)，并集成得到全局坐標(biāo)系下精細(xì)化子結(jié)構(gòu)的邊界位移轉(zhuǎn)換關(guān)系：

式中，C為全局坐標(biāo)系下子結(jié)構(gòu)的邊界轉(zhuǎn)換矩陣：

式中，m1和m2分別為第1 和第2 個(gè)邊界的節(jié)點(diǎn)數(shù)。

根據(jù)主、子結(jié)構(gòu)邊界反力所作虛功相等，圖3(d)中邊界反力之間滿足：

將式(20)代入式(22)，得到屈服構(gòu)件內(nèi)力re：

再將式(23)代入式(3)，得到屈服構(gòu)件的非線性修正力=keue-CTqb。

隔離子結(jié)構(gòu)的增量平衡方程可以表示為：

式中：下標(biāo)i和b分別為隔離子結(jié)構(gòu)內(nèi)部和邊界自由度；ks和q分別為隔離子結(jié)構(gòu)的剛度和外力；Δ為增量。

將式(20)代入式(24)，并將式(24)第二個(gè)等式左乘CT，再通過(guò)靜力凝聚，則有：

4 分析流程

采用數(shù)值子結(jié)構(gòu)方法進(jìn)行結(jié)構(gòu)地震彈塑性時(shí)程分析時(shí)，具體分析流程如下：

1)初始化。令n=0、Fn=0 、=0 、=0、Bs,n=0 以及=0。根據(jù)式(8)集成等效剛度矩陣K，計(jì)算并保存其逆矩陣K-1。根據(jù)總時(shí)步NS P，得到荷載增量 dF=F/NS P。

2)令n=n+1、i=0、=Bs,n-1以 及=，計(jì)算當(dāng)前時(shí)步外力Fn=Fn-1+dF，并由式(9)計(jì)算得到等效外力Pn。根據(jù)式(10)計(jì)算得到主結(jié)構(gòu)位移，更新整體結(jié)構(gòu)屈服構(gòu)件單元狀態(tài)。若屈服單元數(shù)為零，則系統(tǒng)收斂，跳到第7)步；否則，更新布爾矩陣，由式(11)得到塑性自由度位移，進(jìn)入下一步。

5)根據(jù)下式判斷塑性自由度系統(tǒng)的收斂性。

5 算例分析

為了驗(yàn)證本文提出的改進(jìn)的降階牛頓迭代數(shù)值子結(jié)構(gòu)方法(下文簡(jiǎn)稱本文方法)的計(jì)算效率和模擬精度，將其與傳統(tǒng)的牛頓算法(下文簡(jiǎn)稱傳統(tǒng)方法)進(jìn)行對(duì)比分析，同樣以上述15 層3 跨的平面鋼框架結(jié)構(gòu)為例。為避免平臺(tái)的差異影響對(duì)比分析的客觀性，本文算例均在OpenSees 平臺(tái)中進(jìn)行，采用BandGeneral 求解器，計(jì)算硬件條件是Interl Xeon W-2155CPU，主頻是3.3 GHz，內(nèi)存為64 GB。

5.1 屈服單元隔離非線性分析

進(jìn)行動(dòng)力時(shí)程分析之前，將豎向重力荷載分10 步逐漸加到結(jié)構(gòu)上。圖4 給出8 度和9 度罕遇地震作用下整體結(jié)構(gòu)中隔離單元占比時(shí)程曲線。從圖4 中可以看出，結(jié)構(gòu)最終屈服單元數(shù)占比分別是4.4%和9.6%；8 度罕遇地震作用下，1.62 s時(shí)刻結(jié)構(gòu)進(jìn)入彈塑性狀態(tài)，直到12 s 左右屈服單元數(shù)達(dá)到最大；9 度罕遇地震作用下，1.46 s 時(shí)刻進(jìn)入屈服狀態(tài)，直到12 s 左右塑性區(qū)不再擴(kuò)展。

圖4 隔離單元所占比例時(shí)程曲線Fig. 4 Isolated element percentage time-history curves

圖5 給出該平面鋼結(jié)構(gòu)在8 度和9 度罕遇地震作用下的塑性區(qū)發(fā)展歷程及大小，圖中數(shù)值1～6 表示有1 個(gè)～6 個(gè)屈服單元。由圖5 可知，8 度罕遇地震下，1 層～4 層梁端先進(jìn)入屈服狀態(tài)，隨后，5 層～9 層梁端出現(xiàn)塑性區(qū)，其中，1 層～4 層梁端塑性區(qū)面積較大；9 度罕遇地震下，1 層～6 層梁端先進(jìn)入屈服狀態(tài)，然后，7 層～11 層梁端和底層柱端出現(xiàn)塑性區(qū)，其中，1 層～5 層梁端塑性區(qū)擴(kuò)展面積較大，底層柱端塑性區(qū)面積相對(duì)較小。

圖5 塑性區(qū)發(fā)展歷程及大小Fig. 5 Development history and size of plastic zone

5.2 準(zhǔn)確性和可靠性分析

作為對(duì)比，采用傳統(tǒng)方法分析時(shí)，該結(jié)構(gòu)梁柱構(gòu)件同樣劃分30 個(gè)非線性纖維梁?jiǎn)卧鹘y(tǒng)方法和本文方法得到的結(jié)構(gòu)頂點(diǎn)位移和基底剪力時(shí)程曲線如圖6～圖9 所示，可以看出，兩種方法的計(jì)算結(jié)果幾乎一致。

圖6 8 度罕遇地震鋼結(jié)構(gòu)頂層位移時(shí)程曲線Fig. 6 Time-history curves of top displacement for the steel structure under 8 degree rare earthquakes

圖7 9 度罕遇地震鋼結(jié)構(gòu)頂層位移時(shí)程曲線Fig. 7 Time-history curves of top displacement for the steel structure under 9 degree rare earthquakes

圖8 8 度罕遇地震鋼結(jié)構(gòu)基底剪力時(shí)程曲線Fig. 8 Time-history curves of base reaction for the steelstructure under 8 degree rare earthquakes

圖9 9 度罕遇地震鋼結(jié)構(gòu)基底剪力時(shí)程曲線Fig. 9 Time-history curves of base reaction for the steelstructure under 9 degree rare earthquakes

進(jìn)一步，表2 對(duì)比兩種方法的最大頂點(diǎn)位移和最大基底剪力以及它們的相對(duì)誤差。由表2 可知，8 度和9 度罕遇地震作用下，結(jié)構(gòu)最大頂層位移和基底剪力的相對(duì)誤差均不超過(guò)0.2%。由式(28)計(jì)算得到頂層位移時(shí)程曲線的均方差分別為0.14%和0.04%，基底剪力時(shí)程曲線的均方差分別為0.02%和0.01%。進(jìn)一步驗(yàn)證了本文方法的準(zhǔn)確性、可靠性。

表2 傳統(tǒng)方法與本文方法的誤差Table 2 Errors using the traditional and proposed methods

式中：N為總數(shù)據(jù)量；Xi和Yi為第i時(shí)刻響應(yīng)值。

5.3 計(jì)算效率分析

8 度和9 度罕遇地震作用下，傳統(tǒng)方法和本文方法的迭代次數(shù)時(shí)程曲線分別如圖10 和圖11 所示。此外，傳統(tǒng)方法、本文方法以及不動(dòng)點(diǎn)迭代數(shù)值子結(jié)構(gòu)方法(簡(jiǎn)稱原方法)的計(jì)算信息如表3所示。由圖10 和圖11 可知，本文方法的迭代收斂次數(shù)略高于傳統(tǒng)方法，當(dāng)?shù)卣鸷奢d水平增大或迭代收斂誤差減小，本文方法的迭代收斂次數(shù)變化較小。由表3 可以看出，當(dāng)收斂誤差為1.0×10-5mm時(shí)，8 度和9 度罕遇地震作用下，傳統(tǒng)方法最多分別需要4 步和5 步收斂，平均迭代次數(shù)分別為2.7和3.1，本文方法最多分別需要6 步和7 步收斂，平均迭代次數(shù)分別為3.2 和3.7，而原方法最多分別需要24 步和26 步收斂，平均迭代次數(shù)分別為4.6和6.9；當(dāng)收斂誤差為1.0×10-8mm 時(shí)，8 度和9 度罕遇地震作用下，傳統(tǒng)方法最多分別需要5 步和6 步收斂，平均迭代次數(shù)分別為3.2 和3.5，本文方法最多分別需要6 步和7 步收斂，平均迭代次數(shù)分別為3.9 和4.1，而原方法最多分別需要41 步和44 步收斂，平均迭代次數(shù)分別為7.1 和11.1。由此可見(jiàn)，本文方法所需的荷載步最大收斂次數(shù)和平均迭代次數(shù)略大于傳統(tǒng)方法，但明顯小于文獻(xiàn)[10 - 13]中的不動(dòng)點(diǎn)迭代數(shù)值子結(jié)構(gòu)方法。

表3 不同方法迭代信息統(tǒng)計(jì)表Table 3 Table of iteration information using various methods

圖10 8 度罕遇地震迭代步數(shù)時(shí)程曲線Fig. 10 Iteration number time-history curves under an 8-degree rare earthquake

圖11 9 度罕遇地震迭代步數(shù)時(shí)程曲線Fig. 11 Iteration number time-history curves under a 9-degree rare earthquake

當(dāng)收斂誤差為1.0×10-5mm 時(shí)，8 度和9 度罕遇地震下傳統(tǒng)方法和本文方法的計(jì)算信息統(tǒng)計(jì)表如表4 所示?？梢钥闯觯簜鹘y(tǒng)方法分別需要對(duì)整體結(jié)構(gòu)剛度K進(jìn)行7282 次和8267 次集成和分解，計(jì)算耗時(shí)分別為313 min 和398 min，而本文方法只需在初始化過(guò)程中進(jìn)行一次集成和分解，非線性分析過(guò)程中僅需對(duì)Hs矩陣進(jìn)行8409 次和9743 次集成和分解，計(jì)算耗時(shí)分別為157 min和332 min。

表4 傳統(tǒng)方法與本文方法的計(jì)算信息統(tǒng)計(jì)表Table 4 Computation information using the traditional and proposed methods

圖12 給出了傳統(tǒng)方法和本文方法求解矩陣維數(shù)的時(shí)程曲線。從圖12 可以看出：隨著地震荷載水平的增大，塑性自由度增多，迭代更新和分解的Hs矩陣規(guī)模變大；8 度和9 度罕遇地震下，本文方法求解的最大Hs矩陣的維數(shù)分別為666 和1713，而傳統(tǒng)方法求解的K矩陣維數(shù)高達(dá)9315。由此可見(jiàn)，結(jié)構(gòu)非線性程度越低，本文方法求解的Hs矩陣規(guī)模越小，計(jì)算效率越高。

圖12 矩陣運(yùn)算規(guī)模時(shí)程曲線Fig. 12 Matrix operation size time-history curves

基于變剛度的傳統(tǒng)方法需要重復(fù)集成和分解整體結(jié)構(gòu)的剛度矩陣K，以及更新所有構(gòu)件單元的狀態(tài)，而采用本文方法進(jìn)行非線性分析時(shí)，僅對(duì)剛度矩陣K進(jìn)行一次運(yùn)算，迭代過(guò)程中只更新規(guī)模較小的Hs矩陣，且只對(duì)數(shù)量較少的屈服構(gòu)件單元進(jìn)行狀態(tài)更新?？梢酝茰y(cè)，針對(duì)局部非線性的大規(guī)模結(jié)構(gòu)系統(tǒng)的彈塑性分析，本文方法較傳統(tǒng)方法具有顯著的計(jì)算優(yōu)勢(shì)。

6 結(jié)論

本文對(duì)文獻(xiàn)[10 - 13]中不動(dòng)點(diǎn)迭代數(shù)值子結(jié)構(gòu)算法作進(jìn)一步改進(jìn)，提出了降階的牛頓迭代數(shù)值子結(jié)構(gòu)方法，并利用該方法對(duì)一平面15 層3 跨鋼框架結(jié)構(gòu)進(jìn)行地震彈塑性時(shí)程分析，得出以下結(jié)論：

(1)改進(jìn)的降階牛頓迭代數(shù)值子結(jié)構(gòu)方法與傳統(tǒng)牛頓算法得到的結(jié)構(gòu)頂層位移和基底剪力時(shí)程曲線一致，從而驗(yàn)證了本文方法的準(zhǔn)確性和可靠性。

(2)本文方法收斂速率遠(yuǎn)高于不動(dòng)點(diǎn)迭代數(shù)值子結(jié)構(gòu)方法，接近傳統(tǒng)牛頓算法的二次收斂。當(dāng)?shù)卣鸷奢d水平增加或迭代收斂條件變得苛刻，不動(dòng)點(diǎn)迭代數(shù)值子結(jié)構(gòu)方法的收斂次數(shù)明顯增多，而本文方法卻沒(méi)有顯著變化。

(3)本文方法的計(jì)算效率高于傳統(tǒng)牛頓算法。對(duì)于局部非線性結(jié)構(gòu)系統(tǒng)，結(jié)構(gòu)模型塑性自由度規(guī)模較小，本文方法求解的Hs矩陣維度遠(yuǎn)小于傳統(tǒng)方法求解的K矩陣，計(jì)算效率較高，且隨著結(jié)構(gòu)規(guī)模的增加，計(jì)算效率優(yōu)勢(shì)將更加明顯。

(4)數(shù)值子結(jié)構(gòu)方法可以有效提高局部非線性結(jié)構(gòu)系統(tǒng)的模擬精度。隔離屈服構(gòu)件單元并建立精細(xì)化子結(jié)構(gòu)有限元模型，既可以精確捕捉局部屈服構(gòu)件的非線性行為，又可以有效揭示整體結(jié)構(gòu)的失效路徑和塑性分布。