張 東 劉嘯奔 孔天威 楊 悅 武學健 吳 鍇 張 宏

1.中國石油大學(北京)油氣管道輸送安全國家工程研究中心,北京,102249

2.合肥通用機械研究院有限公司國家壓力容器與管道安全工程技術研究中心,合肥,230031

0 引言

長輸管道作為典型的焊接結構,在焊縫區(qū)和熱影響區(qū)難免會出現(xiàn)裂紋、未焊透、未熔合等一系列焊接缺陷[1-2]。近年來我國發(fā)生了一系列管道環(huán)焊縫失效事故,對環(huán)境及人民生命財產造成難以估計的損失[3-4]。自2008年來,以中緬管道、漠大線、陜京三線等為代表的高鋼級管道(X70、X80)在建成試壓和投產運行初期發(fā)生了30余起環(huán)焊縫開裂和泄漏事故[5]。國內外許多學者、機構針對環(huán)焊縫失效問題開展了大量的研究,提出了許多管道環(huán)焊縫安全評價模型[6-7],其中,現(xiàn)有安全評價模型多以環(huán)向全焊縫應力應變本構關系為輸入。但是目前管道環(huán)焊縫的失效一般是由于軸向載荷引起的管道環(huán)向裂紋開裂,這與輸入的環(huán)向材料屬性是不匹配的。吳鍇[8]利用數(shù)字圖像相關(digital image correlation,DIC)技術對X80管道環(huán)焊縫軸向試樣進行了單軸拉伸試驗,測得了試樣頸縮前的應力應變曲線,與環(huán)向全焊縫試樣真應力應變曲線進行了對比,發(fā)現(xiàn)與環(huán)向試樣相比,焊縫金屬軸向強度明顯更小。由此可知,采用環(huán)向全焊縫試樣獲得的應力應變本構關系來表征環(huán)焊縫的材料特性會高估環(huán)焊縫的實際承載能力。為降低管道環(huán)焊縫的失效概率,如何測定焊縫區(qū)材料軸向真實應力應變本構關系成為保障管道安全運行的關鍵一環(huán)。

國內外針對油氣管道焊縫區(qū)材料應力應變本構關系測定問題已開展了大量的研究[9-11],主要包括半經驗解析公式法、小沖桿試驗和DIC技術等方法。在半經驗解析公式方面,CHEN等[12]根據(jù)典型體積元等效能量原理和材料變形區(qū)域的經典空腔假設推導出一種等效能量壓痕模型,基于壓痕的加載部分曲線反向預測了材料單軸應力應變關系,但是該方法的預測結果與拉伸試驗結果的吻合度較差。TU等[13]基于有限元方法提出了一種含缺口試樣真實應力應變曲線轉化為材料等效應力應變曲線的校正解析方法,但是該方法中需要實時記錄試樣的最小橫截面積,對試驗設備的要求較高,且該方法并未通過試驗驗證。在小沖桿試驗方面,CHEN等[14]開展圓片小沖桿試驗獲得了載荷位移曲線,通過等效能量理論推導和有限元數(shù)值模擬相結合的方法提出了一種材料應力應變關系半經驗公式,但是該方法受到小沖桿試驗裝置和試樣尺寸的限制,尚沒有一個統(tǒng)一的半經驗關聯(lián)公式?；谏鲜瞿芰康刃г?劉剛等[15]采用新型長條形小沖桿試樣提出了關于該試樣的Chen-Cai半解析方程的4個系數(shù),進而提出了依據(jù)長條形小沖桿試樣的載荷位移曲線估算金屬材料的單軸本構關系的新型試驗方法,但是該估算方法結果與單軸拉伸試驗結果對比存在不小的誤差。ZHONG等[16]開展小沖桿試驗得到了載荷位移曲線,通過數(shù)據(jù)庫方法確定了Ludwik硬化本構模型的參數(shù),該方法只能計算試樣頸縮前的真實應力應變曲線,且與拉伸試驗結果對比的誤差不可忽略。另外還有學者利用DIC技術探究了非均質材料本構關系。DIC技術通過比較試件變形前后表面數(shù)字圖像差異來實現(xiàn)高精度變形測量,并逐漸應用于焊縫區(qū)域的本構關系測量。REYNOLDS等[17]于1999年發(fā)表了關于使用DIC技術研究焊接區(qū)域本構關系的報告,并基于均勻應力假設對鋁合金板的焊接接頭進行了拉伸試驗,獲得了焊接區(qū)域的板材、熱影響區(qū)和焊縫的力學參數(shù)。但LOCKWOOD等[18]利用有限元的方法對上述實驗結果進行驗證時發(fā)現(xiàn),基于均勻應力假設的結果和有限元模擬存在一定誤差。武旭[19]使用DIC技術并結合MTS拉伸機輸出的載荷獲取了X80管道母材、熱影響區(qū)和焊縫區(qū)材料的應力應變關系。但是DIC技術中試樣真應力是基于體積不變原理得到的,因此在試樣發(fā)生頸縮后,DIC技術得到的應力應變曲線不可靠。

針對以往學者研究的不足,本文基于焊縫區(qū)含缺口圓棒單軸拉伸試驗結果,通過ABAQUS非線性有限元軟件對拉伸試驗進行建模重構,利用貝葉斯正則化反向傳播(Bayesian regularization back-propagation,BRBP)神經網絡和灰狼優(yōu)化算法(grey wolf optimizer,GWO)對ABAQUS中輸入的真實應力應變曲線進行迭代優(yōu)化反演,準確得到了焊縫區(qū)的真實應力應變本構關系。該方法彌補了現(xiàn)有測試高鋼級管道焊縫區(qū)材料應力應變本構關系的不足,以期為高鋼級管道環(huán)焊縫安全評價提供準確的應力應變本構關系,進一步保障油氣管道的安全運行,具有一定的實際工程意義。

1 缺口圓棒拉伸試驗

1.1 試驗設計及流程

為探究高強鋼管道焊縫區(qū)材料真實應力應變本構關系,本節(jié)選取中俄東線現(xiàn)場截取的API-5L X80管線鋼某道焊口作為試驗材料,該道焊口采用全自動焊工藝(熔化極氣體保護實心焊絲電弧焊),未服役且不存在腐蝕等缺陷,試樣取材位置焊材性能合格。以環(huán)焊縫為中間區(qū)域沿管道軸向制作帶缺口圓棒試樣,試樣取樣于環(huán)焊縫填充焊區(qū)域,試驗中所有試樣均為同一位置連續(xù)取樣,圖1所示為試樣與環(huán)焊縫的相對位置關系。參考標準ASTM A370—2017《Standard Test Methods and Definitions for Mechanical Testing of Steel Products》對金屬材料標準拉伸試樣的尺寸設計要求,結合拉伸試驗設備對試樣夾持長度的硬性要求,設計了4組不同缺口類型試驗,每組設置2根試樣以避免試驗結果的偶然性,具體試樣尺寸如圖2所示,其中R為試樣的最大半徑,R1為缺口圓的半徑,R2為試樣的最小半徑,h為缺口圓圓心與試樣表面的距離,試樣夾持段為40 mm,引伸計測量區(qū)長度為20 mm,試樣最大截面直徑為10 mm,其余變化尺寸見表1。

表1 各組試樣尺寸詳細表

圖1 試樣與環(huán)焊縫的相對位置示意圖

圖2 試樣尺寸示意圖

采用GNT200微機控制電子萬能試驗機開展拉伸試驗,選擇合適的夾具,將處理好的試樣裝夾在試驗機液壓夾具上,夾持段取圓棒試樣兩端各40 mm,試驗溫度為室溫20 ℃,引伸計初始長度為20 mm,量程為5 mm。試驗采用位移控制的加載方式,為保證試驗處于準靜態(tài)拉伸狀態(tài),加載速率保持為0.5 mm/min,試驗加載至試樣完全斷裂,試驗設備與加載方式如圖3所示。

(a)試驗設備 (b)試樣加載

1.2 試驗結果

通過焊縫區(qū)含缺口圓棒單軸拉伸試驗得到了4組不同缺口尺寸試樣的載荷-位移曲線。A1-A2試樣、B1-B2試樣的載荷位移曲線分別見圖4a、圖4b,可以發(fā)現(xiàn)A組試驗曲線在位移為0.4 mm之前吻合度較好,在位移為0.4 mm之后A1試樣曲線存在明顯的轉折,該轉折點在宏觀上表現(xiàn)為試樣起裂,A2試樣起裂位移為0.7 mm,A組兩根試樣的起裂位移不同的原因在于焊縫區(qū)域材料的局部差異性。試樣起裂后涉及到材料損傷問題,需要在模型中引入GTN等損傷本構模型,這不在本方法的討論范圍內。B組試驗存在同樣的現(xiàn)象,B1試樣起裂位移為1.0 mm,B2試樣起裂位移為0.7 mm,A組與B組試驗起裂點位移不同的原因在于試樣的最小截面半徑不同。C1-C2試樣、D1-D2試樣的載荷-位移曲線分別見圖4c、圖4d,可以發(fā)現(xiàn)C組、D組試驗曲線在彈性階段基本一致,當試樣局部區(qū)域進入塑性階段后,D1-D2試樣曲線出現(xiàn)了明顯的差異,相較而言,C1-C2試樣曲線差異較小,但仍存在一定的差異性,其原因將在下文進一步討論。

(a)試驗A

觀察4組完成拉伸試驗的圓棒試樣缺口宏觀形貌(圖5)可以發(fā)現(xiàn),A組、B組試樣的變形集中在缺口位置,缺口外區(qū)域沒有發(fā)生明顯的塑性變形;而C組、D組試樣除缺口位置發(fā)生集中變形外,缺口附近的母材區(qū)也發(fā)生了明顯的頸縮現(xiàn)象,且D組試樣頸縮現(xiàn)象較C組更為明顯。由試樣尺寸可知,C組試樣的最小截面直徑為7 mm,D組試樣的最小截面直徑為8 mm,即最小截面面積越大,頸縮現(xiàn)象越明顯。同時,上述試驗宏觀形貌也與各組試驗載荷-位移曲線的差異性結論相一致?？紤]到該缺口圓棒試樣為焊縫與母材組合而成,而該環(huán)焊縫材料處于高強匹配狀態(tài),即焊縫材料屈服強度高于母材的屈服強度,因此,當焊縫區(qū)最小截面與母材最大截面的面積之比大于一定值時,母材區(qū)域將會發(fā)生頸縮現(xiàn)象。并且,隨著最小截面面積的增大,母材的頸縮現(xiàn)象會逐漸明顯,直至頸縮完全出現(xiàn)在母材區(qū)域。

圖5 試樣宏觀形貌圖

2 有限元模型

2.1 邊界與網格

由于無法通過試驗數(shù)據(jù)得到準確的焊縫材料應力應變本構關系,因此,采用非線性有限元軟件ABAQUS建立含缺口圓棒單軸拉伸數(shù)值仿真模型對上述試驗進行重構反演,考慮到模型的對稱性以及后續(xù)優(yōu)化反演算法的效率,該模型采用二維軸對稱模型。模型尺寸與試驗試樣尺寸相同,對缺口右端位置施加YSYMM對稱約束,對模型夾持區(qū)施加Y方向的位移,約束X方向的平動自由度以及UR3方向的轉動自由度。在距離模型右端10 mm處截面設置一個點集,起到試驗中引伸計的作用(圖6)。

圖6 有限元模型示意圖

模型單元采用二次減縮積分單元(CAX8R),對試樣模型進行網格劃分時,考慮到缺口位置應變集中的現(xiàn)象,同時為了提高模型的計算準確度與收斂性,因此在缺口位置與引伸計測量區(qū)域對網格進行加密處理,設定網格單元的尺寸為0.1 mm。采用Y字形過渡網格來減少模型網格數(shù)量(圖6),模型其余部分網格尺寸為1 mm,整個模型包含的單元數(shù)量為3163個。

2.2 材料本構模型

管線鋼材料的應力應變本構模型在工程中存在多種描述方式,在各向同性硬化情況下,通常采用基于冪指數(shù)函數(shù)的硬化模型。NIMA等[20]使用Hollomon模型對管材的本構模型進行表征,發(fā)現(xiàn)冪律模型與管材的真實本構模型較為相近,但該本構模型描述中缺乏彈性模量、屈服強度等管線鋼典型參數(shù)。JIANG等[21]、ZHANG等[22]采用Ramberg-Osgood(RO)模型來描述管線鋼的真應力應變本構模型,該模型將彈性模量、屈服強度納入本構模型中,并通過試驗數(shù)據(jù)與RO模型數(shù)據(jù)對比驗證了本構模型的準確性。

基于上述分析,本文采用油氣管道設計與評價領域常用的RO模型來描述材料的應力應變本構關系,其表達式如下:

(1)

式中,ε為管材應變;σ為管材應力,MPa;E為彈性模量,MPa;α為屈服偏移量;n為硬化指數(shù);σy為管材屈服強度,MPa。

針對母材區(qū)域本構模型,當焊縫區(qū)最小截面與母材區(qū)最大截面的面積之比小于一定臨界值時,母材區(qū)域等效應力始終小于材料屈服強度,即母材區(qū)域總是處于彈性階段。該比例臨界值與材料的匹配系數(shù)有關,匹配系數(shù)越小,該比例臨界值越接近1。由于本文以焊縫區(qū)材料本構關系為研究目標,因此需要保證模型的變形集中在焊縫缺口區(qū)域,防止母材區(qū)域發(fā)生塑性變形,影響后續(xù)反演方法的準確性。綜上所述,本文對母材材料僅賦予純彈性屬性,其中彈性模量取經驗值2×105MPa,泊松比取0.3。

針對焊縫區(qū)域本構模型,由于缺口的存在,該區(qū)域存在應力、應變集中現(xiàn)象,因此該區(qū)域塑性數(shù)據(jù)影響著整個模型結果的準確程度。由式(1)可知該本構模型由E、α、n、σy4個變量決定,其中彈性模量E取經驗值2×105MPa,其余3個變量均屬于未知參量,3個變量的不同組合形成了不同的塑性數(shù)據(jù),因此,確定焊縫區(qū)應力應變本構關系的關鍵問題就是尋找3個變量的最優(yōu)組合。

3 焊縫區(qū)材料本構關系反演方法

3.1 BRBP神經網絡與GWO算法介紹

反向傳播(back propagation,BP)神經網絡是一種按照誤差逆向傳播算法訓練的多層前饋神經網絡。BP神經網絡無論在網絡理論還是在性能方面都已比較成熟,它突出的優(yōu)點就是具有很強的非線性映射能力和柔性的網絡結構。網絡的中間層數(shù)、各層的神經元個數(shù)可根據(jù)具體情況任意設定,并且隨著結構的差異其性能也有所不同[23-25]。但是BP神經網絡也存在學習速度慢、容易陷入局部極小值等缺點。為了克服BP神經網絡的上述問題,貝葉斯正則化BP(BRBP)模型是在BP神經網絡基礎上引入貝葉斯正則化規(guī)則的模型,可提高網絡的訓練速度以及泛化性能,修改后的目標函數(shù)加入了神經網絡權值,其表達式如下:

Ft=λEW+βED

(2)

其中,EW為網絡權重平方和;ED為網絡響應與目標值的殘差平方和;λ、β為目標函數(shù)正則化系數(shù),其大小會直接影響到網絡的擬合效果。若λ?β,則通過調節(jié)網絡規(guī)?？梢缘玫嚼硐氲挠柧氄`差,但可能會出現(xiàn)訓練過度的情況,導致產生過擬合現(xiàn)象;若λ?β,則連接權重可能會大規(guī)模降低,使得網絡規(guī)模減小并簡化,進而出現(xiàn)欠擬合現(xiàn)象,導致訓練誤差達不到預期效果[26-27]。

為了得到最優(yōu)的λ和β,基于貝葉斯分析方法,認為網絡權重和訓練樣本的先驗分布均服從高斯分布。基于后驗函數(shù)最大化的原則,對BR性能函數(shù)進行求解,得到最優(yōu)的正則化系數(shù)λ*和β*,其表達式分別如下:

(3)

(4)

式中,φ為有效權值個數(shù),反映網絡實際規(guī)模;m為神經網絡訓練樣本集的數(shù)目;ω*為神經網絡性能指數(shù)最小時的權值矩陣[28]。

灰狼優(yōu)化算法(GWO)是一種基于灰狼捕食獵物行為開發(fā)的群智能優(yōu)化算法,它具有較強的收斂性能、參數(shù)少、易實現(xiàn)等特點,近年來受到了學者們的廣泛關注,已被成功地應用到了車間調度、參數(shù)優(yōu)化、圖像分類等領域中[29]。整個狼群按照適應度函數(shù)值由大到小依次分為α狼、β狼、δ狼和ω狼四個社會等級。狼群捕獵分為包圍、捕獵、攻擊三個過程。在尋優(yōu)過程中,灰狼個體與獵物之間的距離以及灰狼位置的更新可分別表示為

D=|CXP(t)-X(t)|

(5)

X(t+1)=XP(t)-DA

(6)

其中,D為灰狼與獵物間的距離;t為當前迭代次數(shù);XP、X分別為獵物的位置向量和灰狼的位置向量;A、C為系數(shù)向量,其計算公式分別如下:

A=2μr1-μ

(7)

C=2r2

(8)

其中,μ為收斂因子,隨著迭代次數(shù)從2線性減小到0;r1、r2為[0,1]區(qū)間的隨機向量。經過移動后灰狼群向α狼移動,移動方向由自身位置和隨機向量C決定,移動步長由灰狼與獵物的距離和系數(shù)向量A決定,系數(shù)向量A與C的存在很好地避免GWO算法陷入局部最優(yōu)[30],圖7為灰狼位置更新示意圖。

圖7 灰狼位置更新圖

3.2 反演流程及結果

為探究焊縫區(qū)材料的真實應力應變本構關系,本文基于MATLAB-PYTHON-ABAQUS開展聯(lián)合仿真。以MATLAB為主程序平臺,利用PYTHON對ABAQUS有限元模型進行參數(shù)化計算以及后處理,以Excel為數(shù)據(jù)載體實現(xiàn)三個軟件間的數(shù)據(jù)無損精準傳輸,從而成功地將三個功能強大的計算軟件有機地形成一種智能交互關聯(lián)關系,真正做到一鍵操作即可自動建模、自動提取結果、自動開展優(yōu)化反演分析,顯著提高了研究人員的工作效率,提高了計算結果的精度與準確性,反演流程如圖8所示。

圖8 反演流程圖

基于含缺口圓棒單軸拉伸試驗數(shù)據(jù)與有限元模型開展了焊縫區(qū)應力應變本構關系反演分析。由式(1)可知RO本構模型包含三個可變參數(shù)(σy,α,n),吳鍇[8]利用DIC技術分別對中俄東線X80管道環(huán)焊縫軸向和環(huán)向屈服強度進行了測試,得到材料的屈服強度在630 MPa左右。武旭[19]同樣采用DIC技術對X80管道環(huán)焊縫真應力應變曲線進行了測定,通過擬合方法得到了RO本構模型參數(shù),參數(shù)α的范圍為0.62～0.82,參數(shù)n的范圍為9.2～10.2,屈服強度σy的范圍為542～579 MPa,但是該研究指出其試驗材料為低強匹配材料。考慮到本次拉伸試驗材料為高強匹配材料,結合國家標準GB/T 50470—2017附錄C中各鋼級的參數(shù)取值以及各文獻中得到的RO本構模型參數(shù)取值,保守考慮后基本可以確定參數(shù)α、n的取值范圍。由于筆者前期對該道焊口開展了全焊縫環(huán)向圓棒單軸拉伸試驗,可知該材料的環(huán)向屈服強度約為675 MPa,考慮到環(huán)焊縫材料軸向屈服強度通常小于環(huán)向屈服強度,基于此假定軸向屈服強度范圍為600～700 MPa。其中5組基礎工況參數(shù)取值為基于經驗與文獻的自主取值,具體如表2所示。

表2 RO本構模型初始參數(shù)設置及取值范圍

選擇B1組試驗載荷位移曲線作為目標,將0.05～0.9 mm區(qū)間的曲線均分為24段即可生成25個關鍵點(圖9)。利用非線性有限元軟件ABAQUS分別對5組RO本構模型開展計算得到相應的5組載荷-位移曲線,利用插值方法進一步得到5組載荷-位移曲線的關鍵點信息。以RO本構模型參數(shù)為自變量,以每條載荷-位移曲線25個關鍵點的載荷數(shù)據(jù)信息為因變量,采用BRBP神經網絡進行回歸訓練,訓練集取80%,測試集取20%,允許誤差為10-4,迭代次數(shù)設置為2500,隱含層數(shù)目為20,從而得到RO本構模型參數(shù)與關鍵點的映射函數(shù)關系。

圖9 25個關鍵點選取示意圖

通過BRBP神經網絡擬合映射函數(shù)與拉伸試驗數(shù)據(jù)構建相對均方誤差(relative mean square error,RMSE),其表達式如下:

(9)

其中,Fsim為25個關鍵點處神經網絡輸出的載荷數(shù)據(jù),Fexp為25個關鍵點處試驗曲線中的載荷數(shù)據(jù)。以相對均方誤差為GWO算法的適應度函數(shù)開展優(yōu)化分析,其中,算法灰狼個數(shù)為40,允許誤差為10-4,最大迭代次數(shù)為100。以小于允許誤差和達到最大迭代次數(shù)兩個指標為優(yōu)化算法停止的條件,兩個指標滿足其中一個即跳出優(yōu)化算法。

算法經過18次循環(huán),滿足了算法的允許誤差目標,圖10所示為試驗載荷位移曲線與本文反演方法得到的載荷位移曲線對比,可以發(fā)現(xiàn)二者無明顯的差異,相對誤差僅為0.26%。通過上述優(yōu)化反演方法得到的最優(yōu)應力應變本構曲線及參數(shù)組合如圖11所示。

圖10 B組試驗曲線與反演曲線對比

圖11 最優(yōu)應力應變曲線

分別提取位移d=0.3 mm、0.6 mm、0.9 mm、1.2 mm條件下試樣的等效塑性應變云圖(圖12),可以發(fā)現(xiàn)塑性應變主要集中在試樣缺口處,最大等效塑性應變出現(xiàn)在試樣最小截面的外表面處。由于母材區(qū)域應力在拉伸過程中一直未達到材料屈服強度,因此,在采用本文反演方法過程中,若缺乏母材的材料應力應變本構關系,則可將母材本構關系視為純彈性。但是當試樣最小截面直徑與最大截面直徑之比超過一臨界值時,母材區(qū)域會發(fā)生明顯的塑性變形,因此在未知母材本構關系的情況下利用本方法得到的結果是不可靠的,本文試驗C、試驗D即出現(xiàn)了該現(xiàn)象。

圖12 加載歷程中試樣等效塑性應變分布云圖

3.3 結果驗證

為驗證焊縫區(qū)材料本構關系反演結果的準確性,建立與A組試樣尺寸一致的有限元模型,將基于B組試驗反演得到的最優(yōu)本構模型參數(shù)輸入到有限元模型中,基于母材材料不變開展數(shù)值仿真計算并提取載荷-位移曲線。圖13為A組試驗曲線與有限元曲線的對比圖,可以發(fā)現(xiàn),在彈性階段兩組曲線吻合度較高,在塑性階段兩組曲線存在一定誤差,相對誤差為0.81%,該誤差可能由試樣缺口加工尺寸的微小變化引起,且在工程可接受范圍內。

圖13 A組試驗與有限元曲線對比

4 結論

本文基于MATLAB-PYTHON-ABAQUS聯(lián)合仿真提出了一種通用的管道焊縫區(qū)材料應力應變本構關系優(yōu)化反演方法。開展了4組不同缺口尺寸的拉伸試驗,得到了各試樣的載荷-位移曲線,利用貝葉斯正則化反向傳播(BRBP)神經網絡和灰狼優(yōu)化算法(GWO)得到了焊縫區(qū)材料的真實應力應變本構關系。具體結論如下:

(1)針對焊縫區(qū)含缺口圓棒試樣,在拉伸過程中,需要確保母材區(qū)域不發(fā)生塑性變形,這與試樣的尺寸和材料強度匹配相關。在本文開展的4組試驗中,C組、D組試驗母材區(qū)發(fā)生了明顯的頸縮現(xiàn)象。

(2)以B組試驗曲線為優(yōu)化目標開展優(yōu)化反演,相對誤差為0.26%;將B組得到的最優(yōu)本構模型作為A組有限元模型的輸入,計算得到的載荷-位移曲線與A組試驗得到的曲線的相對誤差為0.81%,充分驗證了反演方法的準確性。

(3)本文提出的思路同樣適用于均質金屬材料大應變范圍應力應變曲線的測定。