一種自緊厚壁圓筒非線性混合硬化模型及殘余應(yīng)力分析

符史仲， 楊國來

(南京理工大學(xué) 機械工程學(xué)院， 江蘇 南京 210094)

0 引言

火炮身管是一種典型厚壁圓筒結(jié)構(gòu)，通過自緊工藝，不僅能提高身管強度和壽命，還能減輕身管質(zhì)量、提高火炮的機動性能。在自緊加工過程中，對殘余應(yīng)力的計算準(zhǔn)確與否將直接影響到自緊身管彈性承載能力及疲勞強度等的設(shè)計計算，因此建立一種能準(zhǔn)確描述身管材料機械特性的本構(gòu)關(guān)系模型至關(guān)重要。

炮鋼是身管的主要材料之一，具有以下材料性能：1)正向加載時應(yīng)力變化較小，近似線彈塑性模型；2)卸載時的彈性模量略低于加載時的彈性模量；3)反向屈服應(yīng)力小于正向屈服應(yīng)力，具有明顯的包辛格效應(yīng)；4)反向加載時應(yīng)力、應(yīng)變具有明顯的非線性關(guān)系。為了描述炮鋼的上述特性，許多學(xué)者提出了相應(yīng)的材料模型。李開平[1]以有限元模擬液壓自緊身管，建立各向異性混合硬化模型，考慮屈服面的平移、縮放和旋轉(zhuǎn)，殘余應(yīng)力計算值與試驗值吻合較好，但該模型的材料參數(shù)較難通過試驗獲取。韓錦章[2]引入包辛格效應(yīng)曲線，提出了基于Prager模型的非線性混合硬化模型，建立自緊厚壁圓筒模型，發(fā)現(xiàn)計算值與試驗值吻合，而由于包辛格效應(yīng)曲線需多次進(jìn)行不同預(yù)應(yīng)變下的拉壓試驗，獲取困難。顏勇等[3]基于Mroz多屈服面模型，建立4段線性模型用于近似描述炮鋼的本構(gòu)關(guān)系模型，分析了自緊厚壁圓筒殘余應(yīng)力分布情況，計算值與試驗值吻合，但多屈服面模型建模復(fù)雜，材料參數(shù)多。馬新謀等[4]推導(dǎo)雙線性隨動硬化數(shù)學(xué)模型，將加載、卸載時的應(yīng)力線性化，描述機械自緊過程的應(yīng)力、應(yīng)變變化，厚壁圓筒殘余應(yīng)力計算值與試驗值較為吻合，然而未考慮的應(yīng)力、應(yīng)變非線性關(guān)系易引起誤差。曾志銀等[5]基于Johnson-Cook模型建立了炮鋼的動態(tài)本構(gòu)模型，并采用高斯- 牛頓迭代法優(yōu)化材料參數(shù)，該本構(gòu)模型驗證試驗采用了Taylor桿試驗系統(tǒng)，試樣長度- 試樣半徑與試驗值最大誤差不超過5%，但未對包辛格效應(yīng)進(jìn)行研究。

由于能模擬材料的平均應(yīng)力松弛和棘輪效應(yīng)，文獻(xiàn)[6]提出的非線性隨動硬化(AF)模型得到廣泛應(yīng)用，且該模型能反映炮鋼的材料特性，可用于身管的有限元模擬中。Chun等[7]基于AF模型和各向異性屈服準(zhǔn)則并增加背應(yīng)力反向加載項，建立了各向異性非線性隨動(ANK)硬化模型，為驗證模型準(zhǔn)確性，進(jìn)行了拉伸仿真，將ANK模型、Chaboche模型、Geng-Wagoner模型等計算得到的應(yīng)力與應(yīng)變曲線與試驗數(shù)據(jù)作比較，結(jié)果表明ANK模型更符合試驗結(jié)果。Kobayashi等[8]建立了循環(huán)黏塑性混合硬化模型，針對臨界面的特征引入了基于背應(yīng)力的分次徑向返回算法，又提出了適合混合硬化模型的新形式一致切線剛度矩陣，并通過缺口棒循環(huán)拉伸仿真驗證了該模型的正確性。謝昌吉[9]建立了超細(xì)晶材料的混合硬化模型，以Al-Mg合金驗證了該模型在不同單調(diào)拉伸和不同循環(huán)拉壓加載工況下的應(yīng)力與應(yīng)變曲線和試驗數(shù)據(jù)吻合。余海燕等[10]基于Barlat屈服準(zhǔn)則建立了混合硬化模型，用U形件的回彈仿真驗證模型正確性，該模型計算值、各向同性硬化仿真結(jié)果均與試驗數(shù)據(jù)進(jìn)行對比，前者回彈量更接近試驗值。

相對于上述模型，自緊身管的殘余應(yīng)力分布受包辛格效應(yīng)的影響較大，且不同的預(yù)應(yīng)變下會出現(xiàn)不同程度的包辛格效應(yīng)[11]，因此在本構(gòu)關(guān)系模型中必須予以考慮。本文建立了一種基于AF模型的非線性混合硬化模型，能描述應(yīng)力、應(yīng)變間的非線性關(guān)系以及隨塑性變形增大而減小的屈服面演變規(guī)律。硬化模型將背應(yīng)力分解成若干項，每個子項的背應(yīng)力均服從AF模型演化規(guī)律，還包含了用于描述包辛格效應(yīng)的非線性各向同性硬化模型。再結(jié)合參量間的本構(gòu)關(guān)系推導(dǎo)出張量形式的一致切線剛度矩陣，以用于整體有限元方程組的迭代求解和提高計算收斂速度。為獲得模型的材料參數(shù)，以拉伸試驗結(jié)果為優(yōu)化目標(biāo)，通過遺傳算法得到最終的材料參數(shù)。以該組參數(shù)為基礎(chǔ)建立自緊厚壁圓筒的有限元模型，分析殘余應(yīng)力分布情況，通過與雙線性硬化模型和試驗數(shù)據(jù)的對比，驗證炮鋼本構(gòu)模型的準(zhǔn)確性。

1 適用于自緊身管的非線性混合硬化模型

針對自緊身管的材料特性，利用有限元軟件Abaqus提供的接口模塊，開發(fā)了一種可適用于自緊身管的非線性混合硬化模型。該模型包含了服從AF硬化模型的背應(yīng)力項和反映屈服面縮放的非線性各向同性硬化項。根據(jù)所建立的本構(gòu)關(guān)系模型，提出了相應(yīng)的本構(gòu)積分算法。為提高整體計算的收斂速度，結(jié)合所建立的本構(gòu)關(guān)系模型和各參量間的彈塑性關(guān)系，推導(dǎo)張量形式的一致切線剛度矩陣。

1.1 非線性混合硬化模型

非線性混合硬化模型包含屈服準(zhǔn)則、流動準(zhǔn)則和硬化準(zhǔn)則。

屈服準(zhǔn)則采用von Mises屈服函數(shù)，流動準(zhǔn)則為關(guān)聯(lián)流動準(zhǔn)則，二者即為基于von Mises準(zhǔn)則的J2流動理論[12]，利用了應(yīng)力張量的偏量部分。

硬化準(zhǔn)則有背應(yīng)力硬化準(zhǔn)則和屈服面硬化準(zhǔn)則。背應(yīng)力演化方程基于AF模型，即

(1)

(2)

卸載時，應(yīng)力狀態(tài)有移入屈服面的趨勢[12]，這時有

(3)

式中：f為屈服函數(shù)；σij為應(yīng)力張量；Dijkl為彈性剛度張量。以(3)式判斷應(yīng)力的加、卸載狀態(tài)。

1.2 隱式積分算法

由于1.1節(jié)所建立的本構(gòu)方程中存在大量的微分方程，需要在數(shù)值計算時對這些微分方程進(jìn)行求解。為保證計算的穩(wěn)定性和精確度，本節(jié)采用隱式積分算法求解本構(gòu)方程，并用于后續(xù)的程序編寫。

假設(shè)第n步為當(dāng)前增量步，各個變量值都已算出，并給定第n+1步的應(yīng)變增量張量Δεij，則可通過積分算法求解出滿足本構(gòu)方程的應(yīng)力張量σij. 以下變量用下標(biāo)n表示第n步的值，用下標(biāo)n+1表示第n+1步的值。

(4)

(5)

1.3 張量形式的一致切線剛度矩陣

為將材料的彈塑性本構(gòu)關(guān)系模型引入有限元分析，不僅需要計算每個積分點的應(yīng)力狀態(tài)，還需要提供與材料的彈塑性本構(gòu)關(guān)系及其積分算法密切相關(guān)的一致切線剛度矩陣，供這一增量步的整體方程組平衡的迭代求解[14]。本文將應(yīng)力方程微分，利用硬化模量、屈服函數(shù)及其他彈塑性方程，推導(dǎo)出一致切線剛度矩陣。該矩陣雖為非對稱張量，但仍可利用對稱張量的矩陣形式方便地編寫成程序，且不存在逆矩陣的求解，從而避免了6×6矩陣的復(fù)雜求逆過程，降低了編程難度，減少了計算時間。

由關(guān)聯(lián)流動準(zhǔn)則可知，塑性因子λ的大小與應(yīng)力增量張量dσij在屈服面外法線上的投影呈正比，即

(6)

式中：H′為硬化模量。進(jìn)而根據(jù)流動準(zhǔn)則，得

(7)

對屈服函數(shù)取微分，并聯(lián)合(1)式和(7)式，得

(8)

式中：f為von Mises屈服函數(shù)。消去(8)式中的公共項(?f/?σmn)dσmn，并聯(lián)合屈服函數(shù)，解(8)式得

(9)

由應(yīng)力與應(yīng)變關(guān)系和(7)式，可得

(10)

解(6)式、(9)式和(10)式得

(11)

為方便之后的方程推導(dǎo)，根據(jù)塑性流動方向nij，可設(shè)如下張量：

(12)

(13)

(14)

(15)

(16)

式中：G為剪切模量。對(13)式進(jìn)行微分，有

(17)

故(17)式可化簡為

(18)

1.4 一致切線剛度矩陣的驗證

表1 牛頓迭代法計算的收斂情況

2 材料參數(shù)選定

目前常采用的材料參數(shù)確定方法有兩類：一類是試驗法，進(jìn)行如拉伸、壓縮、剪切或循環(huán)等試驗，以獲取應(yīng)力與應(yīng)變曲線，然后通過擬合該曲線獲得材料本構(gòu)參數(shù)；另一類是參數(shù)優(yōu)化方法，該方法以材料參數(shù)為優(yōu)化變量，以試驗結(jié)果為優(yōu)化目標(biāo)，通過優(yōu)化方法調(diào)整變量并調(diào)用有限元軟件進(jìn)行計算，直至計算值與試驗值一致。這類方法充分利用了有限元仿真與試驗的優(yōu)點，以簡單、技術(shù)成熟的試驗代替高難度的技術(shù)、不成熟的試驗，具有易實現(xiàn)、方便應(yīng)用的特點[10]。對于本文建立的本構(gòu)模型參數(shù)，彈性模量、泊松比、初始屈服應(yīng)力等參數(shù)，可通過簡單拉伸試驗繪制的應(yīng)力與應(yīng)變曲線獲得；而其他材料參數(shù)，即(1)式和(2)式中的C(i)、γ(i)、Q∞、b，由于難以從試驗中獲取，故需要通過優(yōu)化方法獲得。

參數(shù)優(yōu)化的目的在于尋找一組合適的材料參數(shù)，使得計算結(jié)果與試驗結(jié)果間的誤差最小，即滿足以下目標(biāo)函數(shù)：

(19)

為了高效快速地尋找最優(yōu)解，本文采用遺傳算法，并結(jié)合Isight和Abaqus兩種軟件，用以得到最優(yōu)的材料參數(shù)。相比于傳統(tǒng)優(yōu)化算法，遺傳算法不需要導(dǎo)數(shù)等解析性質(zhì)信息，且具有全局尋優(yōu)能力，比較適應(yīng)于高維、多峰、求解不連續(xù)、具有凹凸性的問題。優(yōu)化流程框圖如圖2所示。

根據(jù)以上方法得到了非線性混合硬化模型的材料本構(gòu)參數(shù)C(i)、γ(i)、Q∞、b的優(yōu)化值，如表2所示。

表2 參數(shù)優(yōu)化值

3 自緊厚壁圓筒殘余應(yīng)力有限元分析

自緊工藝主要有液壓自緊、擠擴自緊和爆炸自緊，其中液壓自緊是在身管內(nèi)壁充滿高壓液體，身管從內(nèi)到外產(chǎn)生塑性或彈性變形。卸壓后根據(jù)原有的應(yīng)力狀態(tài)，壁內(nèi)形成塑性和彈性區(qū)域，產(chǎn)生相應(yīng)的殘余應(yīng)力。

3.1 有限元模型

模型采用內(nèi)徑12.515 mm、外徑26.200 mm的光滑厚壁圓筒，內(nèi)壁承受均勻內(nèi)壓788.4 MPa. 其材料參數(shù)為拉伸彈性模量E=206.0 GPa，卸載彈性模量E′=198.0 GPa，泊松比ν=0.27，拉伸屈服應(yīng)力σ0=935.4 MPa.

本文通過Abaqus/Standard程序建立厚壁圓筒的有限元模型，并應(yīng)用子程序UMAT建立炮鋼的本構(gòu)模型，通過設(shè)置*Depvar和*User Material將狀態(tài)變量和材料參數(shù)輸入至Abaqus軟件中進(jìn)行迭代求解。

3.2 數(shù)值計算結(jié)果及分析

數(shù)值計算得到的切向、徑向殘余應(yīng)力與雙線性混合硬化(BLH，即正反向加載時的應(yīng)力與應(yīng)變關(guān)系簡化為線性)模型得到的結(jié)果如圖3所示。

試驗結(jié)果數(shù)據(jù)源于文獻(xiàn)[1]，該試驗采用膛削法測量厚壁圓筒在不同半徑處(共17處)的殘余應(yīng)力，其中厚壁圓筒的尺寸參數(shù)與本文算例一致。由圖3可以看出，在內(nèi)壁處的切向應(yīng)力，本模型計算結(jié)果與試驗結(jié)果相吻合，在其他應(yīng)力區(qū)的計算結(jié)果與試驗數(shù)據(jù)間的誤差較小，基本吻合。在內(nèi)壁處，BLH模型的計算結(jié)果偏大，在其他應(yīng)力區(qū)的計算結(jié)果與試驗結(jié)果相差較大。

4 結(jié)論

為了描述自緊身管的殘余應(yīng)力分布情況，本文建立了一種基于AF模型的非線性混合硬化模型，并以拉伸圓棒和厚壁圓筒有限元模型為算例進(jìn)行驗證。

1)建立非線性混合硬化模型。硬化模型包含若干背應(yīng)力子項，每個子項的背應(yīng)力均服從AF模型演化規(guī)律；還包含了用于描述包辛格效應(yīng)的非線性各向同性硬化模型。通過本構(gòu)模型間的關(guān)系，推導(dǎo)了張量形式的一致切線剛度矩陣，并模擬簡單拉伸試驗驗證了該矩陣的正確性和有效性。

2)計算殘余應(yīng)力分布情況。以簡單拉伸試驗結(jié)果為優(yōu)化目標(biāo)，通過模型參數(shù)優(yōu)化驗證并確認(rèn)了炮鋼的本構(gòu)模型參數(shù)。由厚壁圓筒的數(shù)值計算結(jié)果可知，本文建立的本構(gòu)關(guān)系模型能較準(zhǔn)確地反映殘余應(yīng)力的分布情況。

由于目前該模型僅考慮了準(zhǔn)靜態(tài)條件下的力學(xué)特性，而高溫動態(tài)條件下(如火炮發(fā)射)材料的許多性能都會發(fā)生改變，發(fā)生諸如高溫蠕變、屈服應(yīng)力增大等現(xiàn)象，這需要更進(jìn)一步地完善本構(gòu)關(guān)系模型。