基于子結(jié)構(gòu)的Woodbury非線性分析方法

蘇 璞，李 鋼，余丁浩

(大連理工大學(xué)海岸和近海工程國家重點(diǎn)實(shí)驗(yàn)室，遼寧，大連 116024)

隨著工程結(jié)構(gòu)規(guī)模和復(fù)雜程度的日益增加，非線性分析已經(jīng)成為研究結(jié)構(gòu)損傷演化規(guī)律的重要手段。材料非線性是土木工程領(lǐng)域常見的一種非線性行為，準(zhǔn)確模擬結(jié)構(gòu)材料非線性行為對(duì)結(jié)構(gòu)性能評(píng)估意義重大。

對(duì)于材料非線性問題，數(shù)值求解過程需重復(fù)合成和分解整體剛度矩陣，尤其是自由度規(guī)模比較大時(shí)，計(jì)算成本大幅增加，耗費(fèi)時(shí)間通常難以接受，盡管計(jì)算機(jī)性能的不斷提高在一定程度上緩解了這一問題，但高效的非線性數(shù)值算法依然是眾多學(xué)者追求的目標(biāo)，許多分析方法不斷涌現(xiàn)，如多尺度方法[1－3]，子結(jié)構(gòu)方法[4－6]等。其中，子結(jié)構(gòu)概念最初用于求解線彈性問題[7]，Clough等[8]將其進(jìn)一步擴(kuò)展到非線性問題中，其核心思想是將整體結(jié)構(gòu)劃分為線彈性子結(jié)構(gòu)和非線性子結(jié)構(gòu)，通過靜力凝聚的方式消去彈性子結(jié)構(gòu)的內(nèi)部自由度，大幅降低了系統(tǒng)自由度數(shù)目，提高了計(jì)算效率，之后該方法被眾多學(xué)者進(jìn)一步完善，并應(yīng)用于不同結(jié)構(gòu)的非線性分析[9－12]，但該方法需要提前預(yù)測(cè)非線性發(fā)生的區(qū)域，其發(fā)展受到一定限制，為克服其局限性，Han和Abel[13]提出了自適應(yīng)子結(jié)構(gòu)方法，該方法可根據(jù)結(jié)構(gòu)非線性狀態(tài)動(dòng)態(tài)的劃分彈性和非線性子結(jié)構(gòu)，無需預(yù)判非線性發(fā)生的區(qū)域；Sheu等[14－15]基于多級(jí)剖分策略提出了可考慮非線性區(qū)域變化的子結(jié)構(gòu)分析方法。此外，有學(xué)者基于混合變分原理提出了雙重子結(jié)構(gòu)方法[16－17]，該方法以界面力為未知量，其界面方程規(guī)模一般小于或等于傳統(tǒng)的子結(jié)構(gòu)方法[18]。盡管子結(jié)構(gòu)方法降低了結(jié)構(gòu)整體剛度矩陣的階數(shù)，但在材料非線性分析過程中仍不可避免的需要對(duì)子結(jié)構(gòu)剛度矩陣進(jìn)行分解運(yùn)算，該過程通常耗時(shí)較多，仍會(huì)在一定程度上降低整體效率。

對(duì)于大多數(shù)工程結(jié)構(gòu)，材料非線性行為一般僅發(fā)生在局部區(qū)域，即結(jié)構(gòu)整體剛度矩陣只有少量元素發(fā)生變化，傳統(tǒng)變剛度法須對(duì)整體剛度矩陣實(shí)時(shí)更新和分解，無法充分考慮結(jié)構(gòu)的局部非線性特征。而在數(shù)學(xué)領(lǐng)域，對(duì)于僅有少量元素變化的矩陣，可通過Woodbury公式快速求逆，目前已有許多學(xué)者將該公式應(yīng)用于局部非線性問題[19－22]。Woodbury公式的使用避免了整體剛度矩陣的實(shí)時(shí)更新和分解，僅需對(duì)小規(guī)模的Schur補(bǔ)矩陣進(jìn)行分解，進(jìn)而有效提升計(jì)算效率。為了將Woodbury公式應(yīng)用于一般有限元問題，李鋼等[23]基于材料應(yīng)變分解的思想提出了隔離非線性有限元法，在該方法中，結(jié)構(gòu)的非線性特征由額外設(shè)置的非線性自由度來體現(xiàn)，而整體剛度矩陣始終保持不變，當(dāng)結(jié)構(gòu)非線性規(guī)模較小時(shí)，結(jié)構(gòu)非線性自由度數(shù)也較少，結(jié)合Woodbury公式僅需對(duì)與非線性自由度相關(guān)的小規(guī)模的Schur補(bǔ)矩陣進(jìn)行分解，從而實(shí)現(xiàn)非線性問題的高效求解[24]。然而，Schur補(bǔ)矩陣通常為滿陣，當(dāng)結(jié)構(gòu)非線性規(guī)模比較大時(shí)，Schur補(bǔ)矩陣階數(shù)也隨之增大，對(duì)其進(jìn)行分解的計(jì)算成本大幅增加，導(dǎo)致Woodbury公式的高效性受到限制。

本文以子結(jié)構(gòu)和隔離非線性有限元法為基礎(chǔ)，提出基于子結(jié)構(gòu)的Woodbury非線性分析方法；該方法將整體結(jié)構(gòu)劃分為若干個(gè)子結(jié)構(gòu)的同時(shí)，也將Schur補(bǔ)矩陣分解為若干個(gè)子矩陣(即子結(jié)構(gòu)的Schur補(bǔ)矩陣)，有效降低了Schur補(bǔ)矩陣的規(guī)模；且整個(gè)分析過程僅需對(duì)子結(jié)構(gòu)的Schur補(bǔ)矩陣進(jìn)行分解，改進(jìn)了 Schur補(bǔ)矩陣階數(shù)較大時(shí) Woodbury公式計(jì)算效率降低的不足，進(jìn)一步拓展了Woodbury公式的適用范圍；依據(jù)時(shí)間復(fù)雜度理論建立效率分析模型，并通過數(shù)值算例證明了本文方法的準(zhǔn)確性和高效性。

1 基本理論

1.1 子結(jié)構(gòu)法

子結(jié)構(gòu)法是將整體結(jié)構(gòu)的求解域劃分為若干個(gè)子域，各子域?qū)?yīng)的部分被定義為子結(jié)構(gòu)，子結(jié)構(gòu)間相鄰的部分稱為界面；同時(shí)，整體結(jié)構(gòu)位移自由度隨之分為內(nèi)部自由度和界面自由度，通過靜力凝聚消去內(nèi)部自由度，即可得到一個(gè)僅與界面自由度相關(guān)的控制方程。因此，子結(jié)構(gòu)法是通過“化整為零”的方式縮減問題規(guī)模，其求解流程為：將整體結(jié)構(gòu)劃分為若干個(gè)子結(jié)構(gòu)；裝配子結(jié)構(gòu)得到界面方程；利用界面信息得到每個(gè)子結(jié)構(gòu)的響應(yīng)[8]。

假定某結(jié)構(gòu)求解域?yàn)?，邊界為Γ，如圖1(a)所示；依據(jù)子結(jié)構(gòu)法將整體結(jié)構(gòu)劃分為若干個(gè)子結(jié)構(gòu)，這里以兩個(gè)子結(jié)構(gòu)為例，如圖1(b)所示；其中：子結(jié)構(gòu)1和子結(jié)構(gòu)2對(duì)應(yīng)的求解域分別為?1和?2；Γ1和Γ2分別為?1和?2所獨(dú)有的邊界，而Γ12為?1和?2的公共界面。

圖1 子結(jié)構(gòu)方法示意Fig.1 Diagram of substructuring method

這里將由子結(jié)構(gòu)組成的系統(tǒng)定義為子結(jié)構(gòu)系統(tǒng)[25]?？紤]界面處子結(jié)構(gòu)間的相互作用力以及位移連續(xù)性條件，則子結(jié)構(gòu)系統(tǒng)的平衡方程和位移協(xié)調(diào)方程分別為[17]：

式中：s為子結(jié)構(gòu)編號(hào)；NS為子結(jié)構(gòu)個(gè)數(shù)，這里NS=2；K(s)、u(s)和F(s)分別為子結(jié)構(gòu)s的剛度矩陣、位移向量和荷載向量；λ為拉格朗日乘子；B(s)是與位移自由度相關(guān)的矩陣。由式(1)解得u(s)，并將其代入式(2)，可得到子結(jié)構(gòu)系統(tǒng)的界面方程為：

其中：為子結(jié)構(gòu)s剛度矩陣K(s)的廣義逆[16]，當(dāng)K(s)非奇異時(shí)，。式(3)通常采用共軛梯度法來求解，將λ回代到式(1)，即可得到各子結(jié)構(gòu)的位移響應(yīng)。在整個(gè)求解過程中，僅需對(duì)各子結(jié)構(gòu)剛度矩陣K(s)進(jìn)行更新和分解，避免了對(duì)結(jié)構(gòu)整體剛度矩陣的相應(yīng)運(yùn)算，可有效降低計(jì)算成本。

1.2 基于Woodbury公式的結(jié)構(gòu)分析方法

Woodbury公式(詳見附錄)現(xiàn)已被廣泛應(yīng)用于高效求解結(jié)構(gòu)非線性問題[19－22]，該類方法可統(tǒng)一稱為基于 Woodbury公式的結(jié)構(gòu)分析方法(以下簡(jiǎn)稱Woodbury法)。隔離非線性有限元法作為一種高效的Woodbury法，其以有限元基本理論為基礎(chǔ)，具有通用性好，適應(yīng)范圍廣的特點(diǎn)。

隔離非線性有限元法的核心思想是將材料應(yīng)變分解為線性和非線性兩部分。以圖2所示的單軸應(yīng)力-應(yīng)變關(guān)系作簡(jiǎn)要介紹，假設(shè)某個(gè)增量步內(nèi)，材料非線性狀態(tài)由A點(diǎn)到B點(diǎn)，應(yīng)變?cè)隽喀う鸥鶕?jù)初始初始彈性模量Ee可分解為：

式中：Δε?為線彈性應(yīng)變?cè)隽?；Δ?為非線性應(yīng)變?cè)隽?。與應(yīng)變?cè)隽喀う艑?duì)應(yīng)的應(yīng)力增量Δσ可表示為：

當(dāng)考慮多維材料本構(gòu)關(guān)系時(shí)，上述應(yīng)變分解過程同樣適用，此時(shí)式(4)和式(5)可相應(yīng)寫為：

式中：Δε、Δε?和 Δε?分別為應(yīng)變向量、線彈性應(yīng)變向量和非線性應(yīng)變向量；Δσ為應(yīng)力向量；De為材料初始彈性本構(gòu)矩陣。

圖2 應(yīng)變非線性分解Fig.2 Nonlinear decomposition of material strain

基于上述應(yīng)變分解思想，結(jié)合虛功原理可得到隔離非線性有限元法的控制方程為[23]：

式中：Ke為結(jié)構(gòu)的初始彈性剛度矩陣；Kr為系數(shù)矩陣，代表結(jié)構(gòu)非線性自由度與節(jié)點(diǎn)力之間的關(guān)系；Krr為分塊對(duì)角矩陣，體現(xiàn)結(jié)構(gòu)材料非線性變形信息；Δu和分別為結(jié)構(gòu)的結(jié)點(diǎn)位移向量和非線性應(yīng)變向量；ΔF為結(jié)構(gòu)結(jié)點(diǎn)荷載向量。采用Woodbury公式求解上述控制方程，可得：

2 基于子結(jié)構(gòu)的Woodbury法

2.1 控制方程

依據(jù)有限元法對(duì)圖1(a)所示結(jié)構(gòu)的求解域進(jìn)行離散化，如圖3(a)所示；采用子結(jié)構(gòu)法將整體結(jié)構(gòu)劃分為若干個(gè)子結(jié)構(gòu)，這里以NS=2為例，相應(yīng)子結(jié)構(gòu)系統(tǒng)如圖3(b)所示。圖3中深色部分表示結(jié)構(gòu)在外荷載作用下可能發(fā)生非線性變形的區(qū)域。各子結(jié)構(gòu)結(jié)點(diǎn)分為界面結(jié)點(diǎn)和內(nèi)部結(jié)點(diǎn)，子結(jié)構(gòu)s的結(jié)點(diǎn)位移自由度N(s)相應(yīng)表示為：

圖3 子結(jié)構(gòu)剖分示意Fig.3 Diagram of substructuring partition

將1.2節(jié)中應(yīng)變分解思想應(yīng)用于各子結(jié)構(gòu)，并考慮子結(jié)構(gòu)之間的變形協(xié)調(diào)關(guān)系，可得子結(jié)構(gòu)系統(tǒng)的平衡方程和位移協(xié)調(diào)方程分別為：

式中：

由式(3)、式(15)、式(16)、式(17)和式(18)可知，子結(jié)構(gòu)法在求解界面方程時(shí)，需對(duì)子結(jié)構(gòu)剛度矩陣K(s)進(jìn)行更新的分解，而基于子結(jié)構(gòu)的Woodbury法僅需對(duì)小規(guī)模的矩陣進(jìn)行分解。

2.2 方程求解

本文控制方程的迭代求解方案基于 Newton-Raphson(簡(jiǎn)稱N-R)格式，每一個(gè)N-R迭代步需依次求解式(16)和式(14)，其中式(16)采用共軛梯度法求解，具體迭代格式如圖4所示。

由圖4可知，在共軛梯度法求解Δλ的過程中，由于合成ΔFB和計(jì)算乘積KBpk涉及到子結(jié)構(gòu)Schur補(bǔ)矩陣的分解或回代運(yùn)算，計(jì)算量較大；其余基本是矩陣和向量或向量和向量的乘積運(yùn)算，計(jì)算量相對(duì)較小，可忽略不計(jì)。

求解式(16)得到Δλ后，將其回代到式(14)即可得到各子結(jié)構(gòu)的位移 Δu(s)。對(duì)比式(9)和式(14)可知，傳統(tǒng)Woodbury法需分解整體結(jié)構(gòu)的Schur補(bǔ)矩陣Ksch，而基于子結(jié)構(gòu)的 Woodbury法僅需分解各子結(jié)構(gòu)的 Schur補(bǔ)矩陣(s=1,2,…,NS)。圖5給出了Ksch和的具體形式，由圖5可知，的階數(shù)與整體結(jié)構(gòu)非線性自由度數(shù)相關(guān)，而的階數(shù)僅與子結(jié)構(gòu)s的非線性自由度數(shù)相關(guān)，因此，的階數(shù)一般低于Ksch的階數(shù)，且的階數(shù)隨子結(jié)構(gòu)數(shù)的增多而降低。

圖4 共軛梯度法迭代格式Fig.4 Iterative format of conjugate gradient method

圖5 Schur補(bǔ)矩陣的階數(shù)示意Fig.5 Dimension of Schur complement matrix

此外，通過對(duì)比式(14)、式(15)、式(17)和式(18)可知，回代求解位移 Δu(s)、計(jì)算乘積KBpk和合成ΔFB的計(jì)算過程相差不大。例如對(duì)于乘積KBpk，其計(jì)算過程從右到左依次執(zhí)行，首先，計(jì)算向量pk與矩陣的乘積，即：

圖6 qk 具體計(jì)算步驟Fig.6 Detailed calculation steps of qk

由此，便完成了乘積KBpk的計(jì)算(上述變量均為中間變量)?；卮蠼馕灰痞(s)和合成ΔFB均可分解成上述類似的步驟進(jìn)行運(yùn)算。值得注意的是，由于計(jì)算ΔFB的過程中已經(jīng)對(duì)矩陣進(jìn)行了分解，因此后續(xù)在計(jì)算KBpk和回代求解位移 Δu(s)時(shí)，只需要對(duì)矩陣進(jìn)行回代即可。

3 時(shí)間復(fù)雜度效率模型

時(shí)間復(fù)雜度作為一種僅與算法本身相關(guān)的效率評(píng)價(jià)方法，具體是指算法執(zhí)行過程中所需的計(jì)算工作量，能相對(duì)客觀地評(píng)價(jià)一個(gè)算法的效率[26－27]。

根據(jù)2.2節(jié)論述可知，基于子結(jié)構(gòu)的Woodbury法在求解過程中，合成ΔFB、計(jì)算KBpk和回代求解位移 Δu(s)三項(xiàng)由于涉及到子結(jié)構(gòu) Schur補(bǔ)矩陣的分解或回代運(yùn)算，計(jì)算量較大，其余運(yùn)算計(jì)算成本相對(duì)較小，可忽略不計(jì)。此外，上述三項(xiàng)的計(jì)算過程也相似，以乘積KBpk為例，分別統(tǒng)計(jì)其計(jì)算過程每一步驟的時(shí)間復(fù)雜度并匯總，即可得到計(jì)算KBpk的時(shí)間復(fù)雜度，見表1。合成ΔFB和回代求解位移 Δu(s)的時(shí)間復(fù)雜度也可通過類似的方式得到，本文不再贅述。將上述三項(xiàng)的時(shí)間復(fù)雜度匯總并忽略低階項(xiàng)，即可得到基于子結(jié)構(gòu)的 Woodbury法的時(shí)間復(fù)雜度，見表2。

表1 計(jì)算KBpk的時(shí)間復(fù)雜度Table 1 Time complexity of calculating KBpk

結(jié)合式(22)，則基于子結(jié)構(gòu)的Woodbury法的時(shí)間復(fù)雜度TPro可表示為：

由式(23)可知，界面方程共軛梯度法求解迭代次數(shù)NCG和非線性自由度數(shù)是影響基于子結(jié)構(gòu)的Woodbury法時(shí)間復(fù)雜度的重要參數(shù)，NCG和越大，該方法的時(shí)間復(fù)雜度越高，即效率越低。

4 算例

4.1 工程概況

某7層鋼框架結(jié)構(gòu)，詳細(xì)尺寸如圖7(a)和圖7(b)所示；框架梁柱截面均為工字型，具體信息如圖7(c)所示。每根梁和柱沿長度方向被劃分為5個(gè)單元，共計(jì)2940個(gè)單元，其類型為3結(jié)點(diǎn)Timoshenko梁，結(jié)點(diǎn)自由度總數(shù)為 15456。所有單元截面的纖維數(shù)量均為31根，上、下翼緣各劃分11根纖維，中間腹板劃分9根纖維。在梁、柱節(jié)點(diǎn)處布置集中質(zhì)量，分布情況見表3。鋼材初始彈性模量為Ee=200 GPa，屈服應(yīng)力為σy=300 MPa，應(yīng)力-應(yīng)變關(guān)系為線性硬化模型，硬化模量為Eh=0.02Ee。地震波選取臺(tái)灣集集地震波，輸入方向?yàn)槠矫骐p向，兩個(gè)方向的峰值加速度均調(diào)幅到1.2g。時(shí)間間隔為0.01 s，總分析時(shí)長為30 s。

表2 基于子結(jié)構(gòu)的Woodbury法的時(shí)間復(fù)雜度Table 2 Time complexity of Woodbury approach based on substructuring method

圖7 7層7跨鋼框架Fig.7 Seven-bay seven-story steel frame

表3 梁柱節(jié)點(diǎn)集中質(zhì)量Table 3 Joint masses of structure

4.2 結(jié)果分析

為驗(yàn)證子結(jié)構(gòu)數(shù)對(duì)本文方法計(jì)算結(jié)果的影響，根據(jù)圖8所示的剖分方案將整體結(jié)構(gòu)劃分為2個(gè)子結(jié)構(gòu)，4個(gè)子結(jié)構(gòu)和8個(gè)子結(jié)構(gòu)三種工況，并分別進(jìn)行非線性動(dòng)力時(shí)程分析，最終得到三種工況下的頂層位移時(shí)程曲線如圖9(a)所示。結(jié)果表明：基于子結(jié)構(gòu)的Woodbury法三種工況的計(jì)算結(jié)果基本一致，其相對(duì)誤差在萬分之一以下，即子結(jié)構(gòu)數(shù)對(duì)該方法計(jì)算結(jié)果基本沒有影響。

為進(jìn)一步驗(yàn)證本文方法的準(zhǔn)確性，分別采用Woodbury法和基于子結(jié)構(gòu)的Woodbury法對(duì)該框架結(jié)構(gòu)進(jìn)行非線性動(dòng)力時(shí)程分析，兩種方法計(jì)算得到的頂層加速度和位移時(shí)程曲線分別如圖9(b)和圖9(c)所示。結(jié)果表明：兩種方法的頂層時(shí)程曲線基本重合，說明基于子結(jié)構(gòu)的Woodbury法的計(jì)算精度與Woodbury法基本一致。此外，在單元的非線性狀態(tài)判定過程中，可確定某一單元所有非線性插值點(diǎn)處不為0的非線性應(yīng)變個(gè)數(shù)，此即為該單元的非線性自由度數(shù)；通過統(tǒng)計(jì)某一分析步所有單元的非線性自由度數(shù)之和，即可得到該分析步結(jié)構(gòu)非線性自由度數(shù)，圖10給出了所有分析步結(jié)構(gòu)非線性自由度數(shù)曲線；為直觀反映結(jié)構(gòu)非線性規(guī)模，定義非線性占比γ，其為結(jié)構(gòu)非線性自由度數(shù)和結(jié)點(diǎn)位移自由度數(shù)的比值，即γ=Nr/N。由圖10可知，結(jié)構(gòu)最大非線性自由度數(shù)為1995，此時(shí)非線性占比為γ≈ 1 2.9%，而在整個(gè)分析過程中，結(jié)構(gòu)平均非線性自由度數(shù)為 371，相應(yīng)的平均非線性占比為γ≈ 2.4%。圖11給出了基于子結(jié)構(gòu)的 Woodbury法三種工況下的共軛梯度法迭代次數(shù)曲線，由圖11可知，隨著子結(jié)構(gòu)數(shù)的增加，共軛梯度法迭代次數(shù)也相應(yīng)增加，其中，相比2個(gè)子結(jié)構(gòu)而言，4個(gè)子結(jié)構(gòu)的迭代次數(shù)僅在非線性占比較高的分析步增加明顯，而8個(gè)子結(jié)構(gòu)的迭代次數(shù)則幾乎在每一分析步均有明顯增長。

圖8 子結(jié)構(gòu)剖分方案Fig.8 Partitioning schemes of substructures

圖9 頂層時(shí)程曲線Fig.9 Time-history curves of top floor

圖10 非線性自由度曲線Fig.10 Curves of nonlinear degrees of freedom

圖11 共軛梯度法迭代次數(shù)曲線Fig.11 Curves of iterations of conjugate gradient method

4.3 效率評(píng)價(jià)

考慮剛度矩陣帶寬與結(jié)構(gòu)結(jié)點(diǎn)自由度數(shù)的關(guān)系，Woodbury法一次N-R迭代求解的時(shí)間復(fù)雜度TIS為[27]：

根據(jù)非線性分析過程中每一次 N-R迭代步結(jié)構(gòu)非線性自由度數(shù)(圖(10))以及求解界面方程所需的共軛梯度法迭代次數(shù)(圖(11))，結(jié)合式(24)和式(23)可統(tǒng)計(jì)出 Woodbury法和基于子結(jié)構(gòu)的Woodbury法每一N-R迭代步的時(shí)間復(fù)雜度，表4給出了兩種方法的最大時(shí)間復(fù)雜度和平均時(shí)間復(fù)雜度，圖12給出了基于子結(jié)構(gòu)的Woodbury法三種工況下的時(shí)間復(fù)雜度曲線。

表4 兩種方法的時(shí)間復(fù)雜度對(duì)比Table 4 Time complexity comparison of two methods

圖12 時(shí)間復(fù)雜度曲線Fig.12 Curves of time complexity

由表4可知，基于子結(jié)構(gòu)的Woodbury法的最大時(shí)間復(fù)雜度和平均時(shí)間復(fù)雜度均明顯低于Woodbury法。值得注意的是，隨著劃分子結(jié)構(gòu)數(shù)的增多，基于子結(jié)構(gòu)的Woodbury法的最大時(shí)間復(fù)雜度逐漸降低；而平均時(shí)間復(fù)雜度雖然也有降低的趨勢(shì)，但8個(gè)子結(jié)構(gòu)時(shí)的平均時(shí)間復(fù)雜度要高于4個(gè)子結(jié)構(gòu)，其主要原因可結(jié)合圖10~圖12來說明。通過對(duì)比圖10~圖12可知，基于子結(jié)構(gòu)的Woodbury法的時(shí)間復(fù)雜度隨分析步的變化規(guī)律與非線性自由度數(shù)和共軛梯度法迭代次數(shù)(以下簡(jiǎn)稱迭代次數(shù))類似，這再一次說明非線性自由度數(shù)和迭代次數(shù)是影響本文算法時(shí)間復(fù)雜度的兩個(gè)關(guān)鍵因素；通過進(jìn)一步分析兩個(gè)因素對(duì)每一分析步時(shí)間復(fù)雜度的影響比重可知，在非線性占比較高的分析步，非線性自由度數(shù)起控制作用，此即基于子結(jié)構(gòu)的Woodbury法的最大時(shí)間復(fù)雜度隨著子結(jié)構(gòu)數(shù)的增多而遞減的主要原因；而在非線性占比較低的分析步，迭代次數(shù)起控制作用，此即大多數(shù)分析步中8個(gè)子結(jié)構(gòu)的時(shí)間復(fù)雜度高于4個(gè)子結(jié)構(gòu)的主要原因。

5 結(jié)論

本文基于子結(jié)構(gòu)法有效降低了Schur補(bǔ)矩陣的規(guī)模，同時(shí)將Woodbury公式應(yīng)用于各子結(jié)構(gòu)，僅需對(duì)子結(jié)構(gòu)Schur補(bǔ)矩陣進(jìn)行分解，避免了整體結(jié)構(gòu)Schur補(bǔ)矩陣的相應(yīng)運(yùn)算；依據(jù)時(shí)間復(fù)雜度理論建立效率模型，并通過數(shù)值算例驗(yàn)證了基于子結(jié)構(gòu)的Woodbury法的準(zhǔn)確性和高效性。得到的主要結(jié)論如下：

(1)與 Woodbury法相比，基于子結(jié)構(gòu)的Woodbury法在保證計(jì)算精度的前提下進(jìn)一步提高了Woodbury公式的計(jì)算效率，拓寬了其應(yīng)用范圍；

(2)子結(jié)構(gòu)數(shù)對(duì)基于子結(jié)構(gòu)的 Woodbury法的計(jì)算精度影響極小，但其對(duì)該方法計(jì)算效率影響較大。為保證本文方法具有較好的計(jì)算性能，可根據(jù)“界面結(jié)點(diǎn)自由度數(shù)與結(jié)構(gòu)結(jié)點(diǎn)自由度數(shù)的比值不超過10%”的經(jīng)驗(yàn)準(zhǔn)則來選取合適的子結(jié)構(gòu)數(shù)。

附錄：Woodbury公式

Woodbury公式是用來對(duì)僅有少量元素變化(低秩修正)的矩陣快速求逆的數(shù)學(xué)工具。假設(shè)某個(gè)n×n階的矩陣A可逆，且其逆已知，A經(jīng)低秩修正后得到矩陣B，即B=A+WVT，其中W和V均為n×p階的矩陣，且p遠(yuǎn)小于n。若矩陣A和(I+VTA-1W)可逆，則矩陣B的逆可通過Woodbury公式進(jìn)行高效求解[19]，即：

由于A-1已知，式(A1)的主要計(jì)算量集中于一個(gè)p×p階的矩陣(I+VTA-1W)的合成和求逆，且低秩修正時(shí)，該矩陣的階數(shù)p遠(yuǎn)小于矩陣B的階數(shù)n(即p<