彭大煒, 張世聯(lián)
(上海交通大學(xué),上海200030)
結(jié)構(gòu)極限強(qiáng)度分析的三種有限元解法研究
彭大煒, 張世聯(lián)
(上海交通大學(xué),上海200030)
介紹了弧長(zhǎng)法、阻尼因子法和準(zhǔn)靜態(tài)法等3種計(jì)算結(jié)構(gòu)極限強(qiáng)度的非線(xiàn)性有限元解法及其求解思路。通過(guò)Reckling No.23模型數(shù)值計(jì)算,對(duì)上述3種解法各自的計(jì)算過(guò)程、求解特點(diǎn)和要點(diǎn)進(jìn)行了比較和歸納。為更準(zhǔn)確、高效地運(yùn)用非線(xiàn)性有限元程序計(jì)算分析結(jié)構(gòu)極限強(qiáng)度提出了合理的建議。
弧長(zhǎng)法;阻尼因子法;準(zhǔn)靜態(tài)法;極限承載能力
近年來(lái),結(jié)構(gòu)設(shè)計(jì)的思路已經(jīng)逐步從傳統(tǒng)的許用應(yīng)力設(shè)計(jì)向極限狀態(tài)設(shè)計(jì)發(fā)展[1]。研究結(jié)構(gòu)的極限強(qiáng)度水平,確定其與設(shè)計(jì)載荷水平之間的確切裕度,已成為結(jié)構(gòu)理性設(shè)計(jì)的核心和基礎(chǔ)。
目前,船體結(jié)構(gòu)極限承載能力研究的方法主要有[2]:實(shí)船事故調(diào)查和模型試驗(yàn)法、直接法(解析法)、以Smith法為代表的逐步崩潰法和非線(xiàn)性有限元法。隨著計(jì)算機(jī)技術(shù)的不斷發(fā)展,非線(xiàn)性有限元法現(xiàn)已成為計(jì)算和評(píng)估結(jié)構(gòu)極限承載能力最理想的方法,其分析問(wèn)題的主要流程如圖1所示。
在非線(xiàn)性有限元法中,有3種不同的分析結(jié)構(gòu)極限強(qiáng)度的解法。它們分別是弧長(zhǎng)法、阻尼因子法和準(zhǔn)靜態(tài)法。本文基于ABAQUS大型通用非線(xiàn)性有限元程序,以Reckling-No.23試驗(yàn)?zāi)P蜑閷?duì)象,分別采用上述3種解法進(jìn)行極限強(qiáng)度的數(shù)值分析,并研究不同解法的求解過(guò)程,對(duì)不同解法的應(yīng)用特點(diǎn)進(jìn)行了歸納和總結(jié),且還針對(duì)如何更準(zhǔn)確、高效地運(yùn)用非線(xiàn)性有限元法計(jì)算分析結(jié)構(gòu)極限強(qiáng)度的問(wèn)題提出合理的建議。
從求解本質(zhì)而言,弧長(zhǎng)法和阻尼因子法是以結(jié)構(gòu)非線(xiàn)性靜態(tài)平衡方程式(1)的求解(牛頓拉-普森迭代)為基礎(chǔ)。而準(zhǔn)靜態(tài)方法是以結(jié)構(gòu)非線(xiàn)性運(yùn)動(dòng)方程式(2)的顯式求解(中心差分法)為基礎(chǔ)。
式中:{P}為載荷列陣;{I}為內(nèi)力列陣;[M]為質(zhì)量矩陣;{¨u}為加速度列陣。
以下分別簡(jiǎn)要介紹這3種解法的基本思路。
弧長(zhǎng)法的基本思路是通過(guò)設(shè)置一個(gè)參數(shù)(弧長(zhǎng)l)來(lái)控制平衡方程的增量迭代和收斂。將(1)式寫(xiě)成增量格式:
式中:[KT]為切線(xiàn)剛度矩陣;{Δu}為位移增量;{ΔP}為載荷增量;{R}為殘差力。
在弧長(zhǎng)法中,第i步迭代的載荷增量{ΔP}i,由載荷增量因子Δ λi和參考載荷{Pref}來(lái)控制,即
將(4)式代入(3)式即得到弧長(zhǎng)法的第i步迭代的增量形格式:
圖2 弧長(zhǎng)法示意圖
弧長(zhǎng)法求解時(shí),是以前一步增量計(jì)算得到的平衡點(diǎn)為圓心,以弧長(zhǎng)增量Δli為半徑,通過(guò)牛頓拉-普森迭代搜索下一個(gè)平衡點(diǎn),如圖2所示。每一步的弧長(zhǎng)增量Δli、載荷增量因子Δ λi和位移增量{Δu}i通過(guò)下面的約束方程來(lái)控制:
通過(guò)迭代,直至殘差力在容差{R}i范圍內(nèi)。當(dāng)?shù)趇步迭代完成時(shí),有
由于在弧長(zhǎng)增量Δli中同時(shí)包含了載荷增量{ΔPi}和位移增量{Δui}的信息,運(yùn)用弧長(zhǎng)法能夠全程跟蹤結(jié)構(gòu)在“加載”(Δ λi>0)和“卸載”(Δ λi<0)過(guò)程中的載荷-位移平衡路徑。但在實(shí)際運(yùn)用時(shí),存在跟蹤失敗[3]和收斂性問(wèn)題。目前弧長(zhǎng)法還在不斷地研究改進(jìn)。
阻尼因子法仍采用傳統(tǒng)的載荷步長(zhǎng)增量來(lái)進(jìn)行非線(xiàn)性迭代求解。其思路是通過(guò)在平衡方程式(1)中引入阻尼力項(xiàng)來(lái)求解結(jié)構(gòu)不穩(wěn)定崩潰的問(wèn)題,求解的平衡方程為
其中{Fv}為阻尼力列陣,取決于求解時(shí)結(jié)構(gòu)變形的速度,由廣義節(jié)點(diǎn)速度{v}來(lái)表征:
式中:c為阻尼系數(shù)[4];M*為人工質(zhì)量矩陣;廣義節(jié)點(diǎn)速度{v}={Δu}/Δt。
加載的初始階段,結(jié)構(gòu)尚處于穩(wěn)定狀態(tài),此時(shí)廣義節(jié)點(diǎn)速度{v}很小,故阻尼力項(xiàng){Fv}對(duì)平衡方程式(10)幾乎沒(méi)有影響,可以忽略。隨著載荷的不斷增加,結(jié)構(gòu)趨向于不穩(wěn)定,當(dāng)外載{P}已經(jīng)不能完全由結(jié)構(gòu)內(nèi)力{I}來(lái)平衡時(shí),結(jié)構(gòu)達(dá)到極限狀態(tài),相當(dāng)部分的應(yīng)變能將釋放轉(zhuǎn)化為動(dòng)能,廣義節(jié)點(diǎn)速率{v}迅速增大。此時(shí),阻尼力項(xiàng){Fv}通過(guò)做功消耗釋放的應(yīng)變能,在平衡方程式(10)中起到維持求解系統(tǒng)的“平衡”作用。
由于采用傳統(tǒng)的載荷步長(zhǎng)來(lái)進(jìn)行非線(xiàn)性迭代,阻尼因子法無(wú)法繼續(xù)有效跟蹤結(jié)構(gòu)在“卸載”過(guò)程中的載荷-位移路徑,取而代之的是一條幾乎水平(斜率為0)的直線(xiàn)。通過(guò)考察阻尼力項(xiàng){Fv}為維持系統(tǒng)“平衡”所做的功所占結(jié)構(gòu)應(yīng)變能百分比的歷史變化曲線(xiàn),就能夠確定結(jié)構(gòu)的極限承載能力。相比弧長(zhǎng)法、阻尼因子法的數(shù)值收斂性要好一些,能夠解決的問(wèn)題也更廣泛。
與弧長(zhǎng)法和阻尼因子法的靜態(tài)求解方法不同,準(zhǔn)靜態(tài)法從本質(zhì)上講是一個(gè)結(jié)構(gòu)動(dòng)態(tài)求解的過(guò)程。在求解時(shí),對(duì)式(2)采用中心差分法進(jìn)行顯示的時(shí)間積分,由一個(gè)增量步的動(dòng)力學(xué)條件計(jì)算下一個(gè)增量步的動(dòng)力學(xué)條件,直至求解時(shí)間結(jié)束。準(zhǔn)靜態(tài)法的基本思路就是用慢速加載的動(dòng)態(tài)分析來(lái)模擬靜態(tài)問(wèn)題,所以,求解關(guān)鍵是要設(shè)置合適加載速率。加載速率過(guò)快會(huì)導(dǎo)致求解結(jié)果的局部性(劇烈的結(jié)構(gòu)局部變形),使計(jì)算結(jié)果偏離“準(zhǔn)靜態(tài)”的要求;而過(guò)慢的加載速率意味著較長(zhǎng)的加載時(shí)間,從而使計(jì)算時(shí)間大幅度增加。故分析時(shí),通常會(huì)取從快到慢多個(gè)加載速率進(jìn)行分析比較,以選定一個(gè)合適的加載速率。判斷加載速率是否合適的一個(gè)重要標(biāo)準(zhǔn)就是分析過(guò)程中結(jié)構(gòu)模型的動(dòng)能與其應(yīng)變能之比,一般準(zhǔn)靜態(tài)的要求是小于5%。
由于中心差分法是條件穩(wěn)定的算法[5],在分析時(shí),時(shí)間步長(zhǎng)Δt必須小于穩(wěn)定性限制Δtstable才能保證求解的穩(wěn)定性:
式中:ωmax為結(jié)構(gòu)最高階固有頻率;Δte為結(jié)構(gòu)模型中最小尺寸的單元的穩(wěn)定時(shí)間步長(zhǎng)。Δte與單元的特征尺度Le、彈性模量E和材料密度ρ有關(guān),是穩(wěn)定性限制Δtstable的上界[6]:
在計(jì)算中,由于結(jié)構(gòu)最高階頻率ωmax不易求得,故時(shí)間步長(zhǎng)Δt就取為Δte。
相比弧長(zhǎng)法和阻尼因子法,準(zhǔn)靜態(tài)方法最大的優(yōu)勢(shì)在于采用中心差分法進(jìn)行顯式時(shí)間積分不存在收斂性的問(wèn)題。因此,準(zhǔn)靜態(tài)法能夠很好地求解更復(fù)雜結(jié)構(gòu)崩潰問(wèn)題,如結(jié)構(gòu)的自接觸和材料的失效問(wèn)題。在求解極限狀態(tài)問(wèn)題時(shí),時(shí)間步長(zhǎng)Δt往往較小,準(zhǔn)靜態(tài)“緩慢”加載可能導(dǎo)致求解的機(jī)時(shí)很長(zhǎng),這時(shí)可采用質(zhì)量放大等方法進(jìn)行調(diào)整。
Reckling(1979)[7]采用箱型剖面模擬實(shí)船剖面進(jìn)行了系列總縱極限承載能力試驗(yàn),以研究船體梁在極限狀態(tài)下各種崩潰模式和剖面應(yīng)力分布。本文選取加強(qiáng)筋較多的Reckling No.23號(hào)模型作為分析對(duì)象,模型跨長(zhǎng)l=500 mm。表1列出了模型截面尺寸和材料屬性。
圖3 Reckling有限元模型一階屈曲模態(tài)
表1 Reckling No.23模型截面尺寸和材料屬性
在進(jìn)行極限強(qiáng)度分析時(shí),可將如圖3所示模型的第一階線(xiàn)性屈曲模態(tài)的變形形式作為結(jié)構(gòu)的初始缺陷引入模型網(wǎng)格。
采用弧長(zhǎng)法時(shí),在模型兩端設(shè)置大小相同、方向相反的參考載荷(中垂彎矩),對(duì)應(yīng)各弧長(zhǎng)增量步的載荷與位移的表達(dá)式如式(7)~式(9)所示。得到的端面彎矩-轉(zhuǎn)角曲線(xiàn)如圖4所示:考察曲線(xiàn)峰值點(diǎn)對(duì)應(yīng)的彎矩,即為采用弧長(zhǎng)法得到的極限載荷。
圖4 載荷-位移平衡路徑曲線(xiàn)(弧長(zhǎng)法)
采用阻尼因子法時(shí),在模型端面設(shè)置的載荷應(yīng)超過(guò)結(jié)構(gòu)的極限載荷,一般可取為端面的塑性彎矩。在結(jié)構(gòu)達(dá)到極限狀態(tài)前,采用阻尼因子法和采用弧長(zhǎng)法得到的端面彎矩-轉(zhuǎn)角曲線(xiàn)幾乎重合;當(dāng)結(jié)構(gòu)達(dá)到極限狀態(tài)時(shí),采用阻尼因子法得到的端面彎矩-轉(zhuǎn)角曲線(xiàn)變成一條水平線(xiàn),如圖5(a)所示。圖5(b)、圖5(c)表達(dá)了中垂加載過(guò)程中,阻尼力項(xiàng)所做的功和結(jié)構(gòu)應(yīng)變能隨端部施加的中垂彎矩和端部轉(zhuǎn)角的變化曲線(xiàn)。通過(guò)綜合考察圖5所示的端面彎矩-轉(zhuǎn)角曲線(xiàn)、彎矩-能量曲線(xiàn)和轉(zhuǎn)角-能量曲線(xiàn)來(lái)確定結(jié)構(gòu)的極限承載能力。
圖5 阻尼因子法
采用準(zhǔn)靜態(tài)法時(shí),在模型兩端施加隨時(shí)間光滑變化的強(qiáng)迫位移(轉(zhuǎn)角)。當(dāng)總的強(qiáng)迫位移值一定時(shí),加載速率取決于加載的總時(shí)間。當(dāng)加載時(shí)間為1 s時(shí)的端部強(qiáng)迫位移加載曲線(xiàn)如圖6(a)所示。這里取加載的總時(shí)間為1/100 s、1/20 s和1/1 s由快到慢三個(gè)加載速率進(jìn)行計(jì)算,輸出端面的彎矩反力,得到如圖6(b)所示的彎矩-轉(zhuǎn)角曲線(xiàn)。同時(shí)考察如圖6(c)所示的結(jié)構(gòu)動(dòng)能和應(yīng)變能歷史曲線(xiàn),確定合適的加載速率。
圖6 準(zhǔn)靜態(tài)法
表2列出了三種數(shù)值解法的計(jì)算結(jié)果及其相對(duì)試驗(yàn)值的誤差。計(jì)算表明,雖然三種有限元解法具體的計(jì)算思路和過(guò)程有所不同,但最終得出的結(jié)果基本一致,且與試驗(yàn)值相吻合。下面對(duì)不同解法的應(yīng)用特點(diǎn)與分析的過(guò)程進(jìn)行討論:
表2 數(shù)值計(jì)算結(jié)果與試驗(yàn)值比較
(1)采用弧長(zhǎng)法應(yīng)對(duì)結(jié)構(gòu)模型的規(guī)模加以控制,以避免出收斂性問(wèn)題。
(2)采用阻尼因子法應(yīng)同時(shí)考察載荷-位移曲線(xiàn)和阻尼力項(xiàng)所做的功占應(yīng)變能百分比的變化曲線(xiàn)。當(dāng)結(jié)構(gòu)到達(dá)極限狀態(tài)時(shí),阻尼力項(xiàng)所做的功將迅速增加,甚至超過(guò)結(jié)構(gòu)的應(yīng)變能。這表明阻尼力項(xiàng)在平衡方程中占據(jù)了主導(dǎo)作用,此時(shí)對(duì)應(yīng)的載荷即為結(jié)構(gòu)的極限載荷。
(3)采用準(zhǔn)靜態(tài)方法應(yīng)取由快到慢若干個(gè)加載速率進(jìn)行分析。加載曲線(xiàn)要求其一階和二階導(dǎo)數(shù)都是光滑的,從而避免由于加載不連續(xù)引起的求解波動(dòng)。通過(guò)考察加載過(guò)程中結(jié)構(gòu)動(dòng)能與應(yīng)變能之比來(lái)判斷加載速率是否合適。計(jì)算表明,采用合適的加載速率得到的準(zhǔn)靜態(tài)結(jié)果與靜態(tài)計(jì)算的結(jié)果接近。
本文對(duì)求解結(jié)構(gòu)極限承載能力的非線(xiàn)性有限元法中的3種解法:弧長(zhǎng)法、阻尼因子法和準(zhǔn)靜態(tài)法進(jìn)行了理論介紹,并通過(guò)對(duì)Reckling No.23模型試驗(yàn)的數(shù)值計(jì)算,對(duì)3種解法的求解思路、計(jì)算特點(diǎn)和關(guān)鍵點(diǎn)進(jìn)行了歸納。計(jì)算分析表明,3種解法都是有效的數(shù)值方法,可以相互作為輔助和補(bǔ)充。當(dāng)求解模型規(guī)模較大時(shí),可以考慮采用阻尼因子法和準(zhǔn)靜態(tài)方法進(jìn)行求解;當(dāng)求解涉及復(fù)雜的結(jié)構(gòu)接觸、材料失效等不容易收斂因素時(shí),采用準(zhǔn)靜態(tài)方法可能更為有效。表3總結(jié)了3種非線(xiàn)性有限元解法各自的特點(diǎn)。
表3 3種非線(xiàn)性有限元解法的特點(diǎn)
[1] Jeom Kee Paik,Anil Kumar Thayamballi.Ultimate Limit State Design of Steel-Plated Structures[M].Chichester,U K:Wiley,2003.
[2] ISSC.Report of ISSC Special Task Committee VI.2-Ultimate Hull Girder Strength[R].Proceedings of ISSC 2000[C].Nagasaki,Japan,2000.
[3] 李元齊,沈祖炎.弧長(zhǎng)控制類(lèi)方法使用中若干問(wèn)題的探討與改進(jìn)[J].計(jì)算力學(xué)學(xué)報(bào),1998,15(4):414-422.
[4] ABAQUS Analysis User’s Manual.7.1.1 SOLVING NONLINEAR PROBL EM[S].Using the damping factor.
[5] 王勖成.有限單元法[M].北京:清華大學(xué)出版社,2003.
[6] 莊茁,張帆,岑松.ABAQUS非線(xiàn)性有限元分析與實(shí)例(ABAQUS數(shù)碼工程師系列叢書(shū))[M].北京:科學(xué)出版社,2005.
[7] Reckling K A.Behavior of box girder under bending and shear[R].Proceeding of the International Ship&Offshore Structures Congress(ISSC),Pairs,1979,2:46-49.
Study on Three FEA Methods of Structure Ultimate Strength Analysis
PENG Da-wei, ZHANG Shi-lian
(Shanghai Jiao Tong University,Shanghai 200030,China)
Three algorithms,Arc-Length method,Damping Factor method and Quasi-Static method in nonlinear FEA method of ultimate strength analysis are presented in this paper.The three algorithms have different analysis processes,features and key points,which are discussed and compared with one another through a numeral calculation by Reckling No23 test model.The result of the comparisons gives a reasonable suggestion for better use of nonlinear FEA method of ultimate strength analysis.
Arc-Length method,Damping Factor method,Quasi-Static method,ultimate strength
U661.41
A
圖1 非線(xiàn)性有限元分析流程
2009-11-03;修改稿收到日期:2009-12-25
彭大煒(1984-),男,碩士研究生,主要從事船體結(jié)構(gòu)極限強(qiáng)度分析研究。
1001-4500(2010)02-0001-05