有限元分析軟件ABAQUS在結(jié)構(gòu)彈塑性地震時程分析中的應(yīng)用摘要：ABAQUS軟件是國際上公認的功能強大的通用有限元分析軟件，具有豐富的材料、單元模型庫，能夠處理高度非線性問題，在結(jié)構(gòu)抗震分析領(lǐng)域有著良好的推廣前景。本文介紹了ABAQUS在結(jié)構(gòu)彈塑性地震時程分析應(yīng)用中的兩個問題：（1）傳統(tǒng)上有限元軟件采用的隱式積分算法，求解過程需占用大量的存儲空間，隨著結(jié)構(gòu)自由度數(shù)量和非線性程度的增加，會遭遇到計算瓶頸，而ABAQUS采用的顯式積分算法，能否有效地解決這個問題;（2）ABAQUS中的塑性區(qū)模型能否通過輸入材料的本構(gòu)關(guān)系得到構(gòu)件各時刻的彈塑性變形和內(nèi)力，從而有效地進行結(jié)構(gòu)的彈塑性性能評比。關(guān)鍵詞：ABAQUS彈塑性時程分析隱式算法顯式算法引言ABAQUS是由美國HKS公司研制開發(fā)的國際上公認的功能強大的有限元軟件之一，可以分析各種固體力學、結(jié)構(gòu)力學系統(tǒng)。特別是能夠處理非常復(fù)雜的問題和模擬高度非線性問題。近些年來ABAQUS優(yōu)秀的分析能力和模擬復(fù)雜系統(tǒng)的可靠性使得ABAQUS被各國的工業(yè)和研究機構(gòu)廣泛采用。需進行彈塑性地震時程分析的結(jié)構(gòu)往往自由度比較多，而且在地震作用下，結(jié)構(gòu)構(gòu)件反復(fù)進入彈塑性受力階段，具有高度非線性，傳統(tǒng)上采用隱式算法求解，求解過程需反復(fù)迭代，占用大量存儲空間，極易遭遇計算瓶頸。更為嚴重的是，當結(jié)構(gòu)出現(xiàn)嚴重的剛度突變（塑性鉸模型）時，可能發(fā)生即使細分步長也不能收斂的情況[1]。而ABAQUS采用的顯式算法能有效地解決存儲空間問題，沒有計算瓶頸，且ABAQUS可采用塑性區(qū)模型通過輸入材料的本構(gòu)關(guān)系來進行結(jié)構(gòu)的彈塑性分析，沒有剛度突變問題。所以ABAQUS是目前解決大規(guī)模結(jié)構(gòu)彈塑性時程分析的最佳方案之一。ABAQUS對時程分析問題的處理眾所周知，結(jié)構(gòu)在地震作用下的振動方程式為：（1）式中：[M]、[C]、[K]分別是結(jié)構(gòu)的質(zhì)量矩陣、阻尼矩陣、剛度矩陣。、、分別是結(jié)構(gòu)的加速度向量、速度向量和位移向量，是地面加速度向量。結(jié)構(gòu)地震反應(yīng)時程分析的過程即求解振動方程的過程，ABAQUS采用直接積分法來求解振動方程，ABAQUS/Standard模塊中采用隱式算法進行積分，ABAQUS/Explicit模塊中采用顯式算法進行積分。在彈塑性時程分析中，剛度矩陣是隨時間變化的函數(shù)，選擇合適的彈塑性模型是關(guān)鍵問題，ABAQUS可采用塑性區(qū)模型等。2.1積分算法2.1.1隱式積分算法ABAQUS/Standard應(yīng)用Hilber-hughes方法對運動方程進行隱式的時間積分，Hilber-hughes是Newmark-β法的改進法?；痉匠倘缦拢海?）（3）（4）式中是數(shù)值阻尼系數(shù)，如果=0.0就是Newmark-β法。在ABAQUS中，系數(shù)γ取1/2，β取1/4，所以積分算法無條件穩(wěn)定。采用隱式積分進行靜力分析和彈性動力時程分析時，擬剛度矩陣只需組集一次便可多次求解，求解效率很高。且在幾何非線性和彈塑性非線性分析中，積分步長可以取的比較長，求解速度較快。但擬剛度矩陣求解需消耗很大的存儲空間，其大小隨結(jié)構(gòu)自由度增加呈幾何級數(shù)增長；在求解大規(guī)模問題中，極易遭遇計算瓶頸。在處理高度非線性問題時，需不斷迭代縮小步長求解，每次迭代都需重組剛度矩陣，計算代價很高。2.1.2顯式積分算法ABAQUS/Explicit應(yīng)用中心差分方法對運動方程進行顯式時間積分，動力學平衡方程如下：（5）對時間顯式積分：（6）（7）式中M是質(zhì)量矩陣，Pi是所施加外力，Ii 是單元內(nèi)力。求解過程如下：在當前增量步開始時，用動力學平衡方程(5)計算加速度，對加速度在時間上進行積分時采用中心差分方法，在計算速度的變化時假定加速度為常數(shù)。應(yīng)用這個速度的變化值加上前一個增量步中點的速度來確定當前增量步中點的速度，速度對時間的積分并加上在增量步開始時的位移以確定增量步結(jié)束時的位移。采用顯式積分進行時程分析時，為滿足加速度在時間增量內(nèi)保持不變的假定，積分步長必須小于模型自振周期中的最小值，積分時間步長一般是隱式步長的1/1000～1/100。以加速度為計算的基本要素，在處理實際不產(chǎn)生加速度的靜力加載問題時，需分多步緩慢加載，耗時很長。但其無需進行矩陣求逆，所求解方程均為一元方程，無需迭代求解，每步均可保證收斂。對存儲空間的消耗與單元數(shù)目成正比，資源占用率低，求解速度取決于CPU浮點計算速度，適合處理大規(guī)模問題。隱式算法與顯式算法的比較顯式算法是有條件穩(wěn)定，積分步長依賴于模型的最高階自然頻率，而與載荷類型及持續(xù)時間無關(guān)；隱式算法是無條件穩(wěn)定，積分步長的大小通常參考精度和收斂情況來決定。顯式算法最顯著的特點是無需進行隱式方法中所必需的整體切線剛度矩陣求逆。由于是“顯式”地前推模型狀態(tài)，所以不需要迭代和收斂準則。在隱式分析中，每一次迭代都需要求解大型線性方程組，這一過程需要占用相當數(shù)量的計算資源、磁盤空間和內(nèi)存。顯式算法的求解速度取決于CPU的浮點計算能力，而隱式算法的求解速度取決于虛擬內(nèi)存與磁盤之間的I/O運算。2.2彈塑性模型描述單元（梁、柱）的常用彈塑性模型有塑性鉸模型和塑性區(qū)模型。塑性鉸模型將塑性變形限制發(fā)生在單元的局部區(qū)域，其余部位則永遠保持彈性，這是對工程實際情況的簡化，該模型的單元本構(gòu)關(guān)系大多采用桿件的內(nèi)力-變形關(guān)系，即與桿件截面相關(guān)的本構(gòu)關(guān)系。塑性區(qū)模型則是在單元長度方向和截面內(nèi)逐漸進入塑性，是對工程實際情況的更精確的模擬。該模型的單元本構(gòu)關(guān)系采用材料一點的應(yīng)力-應(yīng)變關(guān)系，單元剛度陣則由材料一點的應(yīng)力－應(yīng)變關(guān)系對截面和單元長度積分得到，因而軸力和雙向彎矩耦合作用的屈服準則和滯回規(guī)律得到較精確的模擬。目前高端的大型通用軟件如ABAQUS、ANSYS等均有塑性區(qū)模型，雖然塑性區(qū)模型計算工作量偏大，但隨著現(xiàn)代微機計算能力的快速提高，使得該模型得以進入實用化階段。ABAQUS有豐富的材料庫，可以提供多種材料本構(gòu)關(guān)系及失效準則模型，如塑性模型有金屬動力硬化模型、混凝土損傷模型等。同時為了方便用戶開發(fā)自己需要的材料模型，ABAQUS提供了方便靈活的二次開放平臺，用戶可以根據(jù)需要編寫用戶子程序（UMAT）。在模擬地震作用下結(jié)構(gòu)的反應(yīng)時，混凝土材料選擇混凝土彈塑性損傷模型，鋼材選擇動力硬化模型（如圖1鋼材雙線型應(yīng)力－應(yīng)變關(guān)系曲線圖1）。算例分析3.1隱式算法與顯式算法的對比一棟10層三跨框架結(jié)構(gòu)[2]，層高均為3.6m，見圖2。各樓層質(zhì)量取為6.21×104kg，梁柱采用Q235鋼材，雙線性動力強化模型，強化系數(shù)0.02，阻尼比0.02。結(jié)構(gòu)構(gòu)件截面參數(shù)見表1。表1構(gòu)件截面參數(shù)構(gòu)件類別寬度/mm高度/mm翼緣厚度/mm腹板厚度/mm中柱5757252016邊柱5006502016所有梁3257002016圖2結(jié)構(gòu)分析模型使用ABAQUS計算得出結(jié)構(gòu)的前三階周期為0.961s、0.389s、0.1705s?？疾旖Y(jié)構(gòu)在地震作用在下的響應(yīng)，選用EI-centro地震記錄。分別采用隱式算法和顯式算法計算20s內(nèi)結(jié)構(gòu)在8度罕遇、9度多遇地震下的結(jié)構(gòu)響應(yīng)。在計算8度罕遇地震結(jié)構(gòu)響應(yīng)時，模型材料采用金屬動力硬化彈塑性模型，在計算9度多遇地震下結(jié)構(gòu)響應(yīng)時，模型材料采用彈性模型。計算所用機器的CPU為P4-3.1G，內(nèi)存1G。用隱式算法計算時，使用ABAQUS/Standard模塊，分別采用B23單元和B21單元，參數(shù)設(shè)置：采用ABAQUS固定積分步長為0.02s。用顯式算法計算時，使用ABAQUS/Explicit模塊，使用B21單元，參數(shù)設(shè)置：積分步長分別采用2e-6，1e-5兩種情況。圖3給出了ABAQUS/Standard和ABAQUS/Explicit得到的8度罕遇地震作用下的底部剪力時程曲線的對比圖（圖中實線表示隱式算法結(jié)果，虛線表示顯式算法結(jié)果）。表2給出的是隱式分析與顯式分析的計算費用對比。(a)隱式算法采用B23單元時(b)隱式算法采用B21單元時圖38度罕遇地震作用下底部剪力時程曲線表2隱式算法與顯式算法的計算費用對比CUP使用率所占內(nèi)存（MB）耗時(s)磁盤空間（MB）隱式算法8度罕遇地震20％～30％820MB197.5s17.2MB9度多遇地震20％～25％750MB193.6s17.2MB顯式算法8度罕遇地震2e-650%~60%670MB5400s1.84MB1e-550%~60%662MB2220s1.77MB9度多遇地震2e-650%~60%665MB5360s1.77MB1e-550%~60%663MB2100s1.80MB從圖3、表2可知：（1）顯式算法中，積分步長越小，由其得到的結(jié)果與由隱式算法得到的結(jié)果越接近。一般地，采用系統(tǒng)默認的自動劃分積分步長足以保證精度。（2）顯式算法采用積分步長2×10e-6時，底部剪力的時程曲線與隱式的在形狀和幅值上基本相同（圖3-a），差別源于在分析中二者分別采用了B21單元和B23單元。實際工程運算時，基于高精度的目的，通常采用B23單元。如果ABAQUS/Standard采用B21單元，那么由兩種算法得到的結(jié)果幾乎沒有差別（圖3-b）。（3）從計算成本上看，隱式算法耗費的存儲空間比較大，而CPU使用率很低。而顯式算法的計算成本較低，有效地使用了CPU資源。隨著自由度數(shù)量的增加，特別是在彈塑性分析時，隱式算法需要大量迭代，方程求解器對資源的需求迅速增加，因此在實際工程運算上，隱式算法的計算瓶頸在于可用存儲空間的大小。（4）由于該模型的單元數(shù)較少，沒有體現(xiàn)出顯式算法在耗時上的優(yōu)越性；在自由度數(shù)目不多，非線性不強烈的情況下，用隱式算法比用顯式算法耗時少。為了比較規(guī)模龐大的結(jié)構(gòu)在顯式算法與隱式算法方面的差異，將上述算例擴展為一個20層15跨的純鋼框架模型，模型其他條件同上例。計算10s內(nèi)8度罕遇地震下的結(jié)構(gòu)響應(yīng)。顯式算法和隱式算法均采用自動劃分積分步長。表3給出的是隱式算法與顯式算法的計算費用對比。表3隱式算法與顯式算法的計算費用對比CPU使用率所占內(nèi)存（MB）耗時磁盤空間（MB）隱式算法20%～30％10243小時15分鐘716顯式算法60％～65％7201小時9分鐘59.2從表3可以看出在結(jié)構(gòu)自由度比較多時，顯式算法比隱式算法耗時明顯減少，存儲空間也小得多，所以在處理結(jié)構(gòu)自由度比較多時，宜采用顯式算法。3.2塑性區(qū)模型與塑性鉸模型的對比某5層辦公樓采用單跨輕鋼框架結(jié)構(gòu)[3]，構(gòu)件截面尺寸見表4，鋼筋混凝土樓面，壓型鋼板加保溫層屋面，墻體為蒸壓輕質(zhì)加氣混凝土(ALC)板，（見圖4）。第2～第4層每個節(jié)點的質(zhì)量為23000kg，頂層節(jié)點質(zhì)量為3410kg。阻尼系數(shù)取0.02，梁柱使用Q235鋼，強化系數(shù)0.02。選用EI-centro地震記錄，時間步長取0.02s，計算抗震設(shè)防烈度分別為8、9度罕遇地震下的結(jié)構(gòu)響應(yīng)。表4構(gòu)件截面參數(shù)構(gòu)件類別寬度/mm高度/mm翼緣厚度/mm腹板厚度/mm柱3005501812樓面梁300600188屋面梁25035086圖4結(jié)構(gòu)分析模型由于模型的自由度較少，使用ABAQUS中的隱式算法進行求解。分別采用B23單元和B21單元進行求解。采用B23單元時每根構(gòu)件分別劃分為1個單元和7個單元兩種情況；采用B21單元時每根構(gòu)件分別劃分為1、2、3、7個單元四種情T1=0.976sT2=0.296sT3=0.157s況。材料采用金屬動力硬化模型（KinematicPlastic）圖5結(jié)構(gòu)的前三個振型，積分步長取0.02s。計算得出結(jié)構(gòu)的前三周期分別為0.976s，0.296s，0.157s。結(jié)構(gòu)的前三個振型如圖5所示。圖6分別給出了由塑性鉸模型[3]和塑性區(qū)模型得到的9度罕遇地震下底部剪力時程曲線，圖7分別給出了由塑性鉸模型[3]和塑性區(qū)模型得到9度罕遇地震下頂點位移時程曲線。8、9度罕遇地震下結(jié)構(gòu)的底部剪力與層間位移曲線如圖8、9所示。圖10給出了由計算結(jié)果得到的進入塑性構(gòu)件的材料本構(gòu)關(guān)系。表5所示的是由塑性區(qū)模型與塑性鉸模型[3]得到的結(jié)果對比。(a)塑性鉸模型(b)塑性區(qū)模型圖69度罕遇地震下的底部剪力時程曲線(a)塑性鉸模型(b)塑性區(qū)模型圖79度罕遇地震下的頂點位移時程曲線.(a)塑性鉸模型(b)塑性區(qū)模型圖88度罕遇地震下的結(jié)構(gòu)底部剪力與層間位移滯回曲線(a)塑性鉸模型(b)塑性區(qū)模型圖99度罕遇地震下的結(jié)構(gòu)底部剪力與層間位移滯回曲線表5塑性區(qū)模型與塑性鉸模型的對比最大位移和剪力8度罕遇9度罕遇塑性鉸模型塑性區(qū)模型差值塑性鉸模型塑性區(qū)模型差值頂點最大位移/ccm19.8815.9219.9%43.7518.557.7%底部最大剪力/kkN656.1377.342.5%727.141043.6%最大層間位移/ccm7.45.328.3%14.636.3556.6%從這些圖表中可看出：（1）塑性區(qū)模型的計算結(jié)果與塑性鉸模型的結(jié)果對比，頂點位移和底部剪力小很多；并且隨著地震最大加速度的提高，位移差距增大，底部剪力差距增加較小，原因可能是因為在塑性區(qū)模型中，在單元長度方向和截面內(nèi)逐漸進入塑性，而在塑性鉸模型中，采用集中塑性鉸假設(shè)，塑性鉸的幾何長度為0。梁端截面進入塑性后即形成塑性鉸，所以在地震反復(fù)運動作用下，塑性鉸模型的塑性變形大些，這可以從圖8和圖9的底部剪力與層間位移曲線中看出來。地震運動越劇烈，在地震作用過程中，塑性變形差距越大，兩個模型間的頂點位移差距越大。而底部剪力增量等于剛度矩陣乘以位移增量，雖然在9度罕遇地震下，塑性鉸模型的位移增量比較大，但由于形成塑性鉸時，剛度比較小，所以底部剪力的差距增加較小。（2）從圖10可以看出，反演的材料應(yīng)力-應(yīng)變關(guān)系與金屬材料的雙線型滯回曲線吻合良好，ABAQUS采用的動力強化塑性模型能很好地模擬鋼材等金屬材料的滯回性能。所以ABAQUS采用塑性區(qū)模型符合工程實際情況。(a)8度罕遇地震(b)9度罕遇地震圖10塑性區(qū)模型反演的材料應(yīng)力-應(yīng)變關(guān)系從圖11可以看出單元對結(jié)果的影響，采用B23單元的模型，網(wǎng)格劃分稀疏對結(jié)果影響不大；采用B21單元的模型，隨著網(wǎng)格劃分加密，時程曲