黃輝，趙耀，袁華，王波

（華中科技大學(xué)船舶與海洋工程學(xué)院，湖北武漢 430074）

0 引言

焊接是現(xiàn)代機(jī)械工業(yè)中的重要連接工藝，是結(jié)構(gòu)產(chǎn)生殘余變形和應(yīng)力的主要來(lái)源。從20世紀(jì)50年代開始，人們開始焊接變形和應(yīng)力的數(shù)值計(jì)算，取得了一系列的研究成果。近10多年來(lái)，隨著計(jì)算機(jī)的迅猛發(fā)展，使焊接過(guò)程的熱彈塑性計(jì)算成為模擬焊接變形、應(yīng)力、組織轉(zhuǎn)變等復(fù)雜現(xiàn)象的主要手段。焊接過(guò)程中材料在高溫條件下發(fā)生分解、相變，而且不同溫度下的材料性能具有很大差異，使焊接局部區(qū)域經(jīng)歷高度非線性的過(guò)程，要模擬整個(gè)結(jié)構(gòu)的焊接變形和殘余應(yīng)力，需要進(jìn)行熱機(jī)耦合非線性有限元分析。由于焊縫區(qū)溫度和應(yīng)力應(yīng)變的梯度很大，進(jìn)行有限元建模分析時(shí)，網(wǎng)格劃分必須足夠精細(xì)，并且進(jìn)行時(shí)間增量步的分析。對(duì)于大型復(fù)雜結(jié)構(gòu)，其焊縫長(zhǎng)度和焊接時(shí)間都很長(zhǎng)，有限元網(wǎng)格規(guī)模和計(jì)算時(shí)間步勢(shì)必急劇膨脹，進(jìn)行完整的熱彈塑性有限元分析將十分困難。

網(wǎng)格重劃分技術(shù)可以處理網(wǎng)格畸變及單元變形過(guò)大的問(wèn)題，使熱軋、巖土等大變形分析得以進(jìn)行下去，并保證計(jì)算結(jié)果的精度和可靠性。S．Brown和H．Song［1］結(jié)合網(wǎng)格重劃分和動(dòng)態(tài)子結(jié)構(gòu)技術(shù)，采用殼體單元分析了平板激光焊接后的變形和應(yīng)力，所使用的時(shí)間與完整的熱彈塑性計(jì)算相比，僅為后者的七分之一。L．－E．Lindgren［2］基于一種可變節(jié)點(diǎn)的體單元，開發(fā)了自適應(yīng)網(wǎng)格技術(shù)，并應(yīng)用三維有限元方法成功模擬了大型銅罐的電子束焊接，大幅縮短了計(jì)算的時(shí)間。本文開展了三維全六面體網(wǎng)格重劃分技術(shù)研究并將其應(yīng)用于焊接數(shù)值模擬，對(duì)平板和圓筒的焊接模型分別進(jìn)行了焊接熱力學(xué)仿真。

1 網(wǎng)格重劃分

焊接過(guò)程中，熱影響區(qū)的溫度、應(yīng)力應(yīng)變梯度較大，具有高度的非線性特征，有限元模擬中單元應(yīng)當(dāng)劃分得很細(xì);然而，遠(yuǎn)離焊縫的絕大部分區(qū)域仍然處于線彈性狀態(tài)，單元可以劃分得相對(duì)較粗。如果不斷對(duì)焊接模型進(jìn)行網(wǎng)格重新劃分，使細(xì)密的網(wǎng)格隨著熱源一起移動(dòng)，就能節(jié)省大量的計(jì)算資源并保證計(jì)算結(jié)果的精度。在網(wǎng)格重劃分前后，必須保證新舊網(wǎng)格的節(jié)點(diǎn)占據(jù)相同幾何空間，以順利地進(jìn)行應(yīng)力應(yīng)變等變量的映射，使分析求解得以連續(xù)進(jìn)行下去。下面介紹節(jié)點(diǎn)遷移和變量映射2個(gè)步驟。

1.1 節(jié)點(diǎn)遷移

由于網(wǎng)格疏密程度在網(wǎng)格重劃分前后變化較大，不得不完全改變節(jié)點(diǎn)和單元的編號(hào)，重新構(gòu)建有限元模型的網(wǎng)格。為了匹配新網(wǎng)格節(jié)點(diǎn)在舊網(wǎng)格中的位置，必須對(duì)網(wǎng)格單元一一尋址，逐一判斷節(jié)點(diǎn)是否處于某單元內(nèi)部［3］。節(jié)點(diǎn)遷移的基本問(wèn)題就是已知單元8個(gè)節(jié)點(diǎn)的坐標(biāo)和位于該單元中1個(gè)節(jié)點(diǎn)的坐標(biāo)，反求這1個(gè)節(jié)點(diǎn)在該單元對(duì)應(yīng)的等參元中的坐標(biāo)。

對(duì)于六面體單元，其形函數(shù)為:

根據(jù)等參元的性質(zhì)，單元中任一點(diǎn)處的整體坐標(biāo)

將表1中坐標(biāo)值代入式（1）和式（2），可得如下方程組:

式中:ai，bi，ci為已知量xi，yi，zi的函數(shù)，它們的值均為常數(shù)。上述高階齊次方程組的求解方法有牛頓法、割線法、布朗方法等，本研究采用了牛頓法，對(duì)于比較規(guī)則的網(wǎng)格單元，進(jìn)行為數(shù)不多的幾次迭代其結(jié)果就能收斂。求得節(jié)點(diǎn)在局部坐標(biāo)系中的坐標(biāo)后，可以用形函數(shù)對(duì)節(jié)點(diǎn)位移進(jìn)行插值，從而得到下一步分析時(shí)的節(jié)點(diǎn)在整體坐標(biāo)系中的坐標(biāo)值。

1.2 變量映射

確定新網(wǎng)格的節(jié)點(diǎn)坐標(biāo)后，需要將舊網(wǎng)格中求解的應(yīng)力應(yīng)變場(chǎng)映射到新網(wǎng)格上，以保證分析的連續(xù)性和準(zhǔn)確性。采用了如下映射算法［4］:

1）采用插值技術(shù)，將上一步網(wǎng)格中積分點(diǎn)上的變量（應(yīng)力、應(yīng)變等）外插到每一個(gè)單元的節(jié)點(diǎn)上，然后對(duì)所有鄰接同一節(jié)點(diǎn)的單元的變量取平均。

2）新網(wǎng)格的積分點(diǎn)的位置由舊網(wǎng)格積分點(diǎn)來(lái)確定，分兩步:①先找到該點(diǎn)所屬舊網(wǎng)格的單元，這樣就可以確定該點(diǎn)在單元中的位置（該步驟假定新網(wǎng)格中所有積分點(diǎn)位于舊網(wǎng)格邊界內(nèi)）。②變量由舊網(wǎng)格的節(jié)點(diǎn)上內(nèi)插到新網(wǎng)格的積分點(diǎn)上。

2 計(jì)算實(shí)例

分析了平板表面堆焊和圓筒環(huán)焊縫的焊接過(guò)程，分別比較了完整模型與網(wǎng)格重劃分模型的計(jì)算結(jié)果。平板堆焊模型的尺寸為500 mm×300 mm×6 mm（圖1）。圓筒的內(nèi)徑為150 mm，外徑為156 mm，長(zhǎng)度為400 mm（圖2）。平板和圓筒的焊接條件見表2。SS400的材料熱物性與熱力學(xué)性能參數(shù)曲線見文獻(xiàn)［5］。

2.1 熱傳導(dǎo)計(jì)算

采用順次耦合的熱彈塑性有限元計(jì)算方法，分別對(duì)平板模型和圓筒模型進(jìn)行焊接熱力學(xué)仿真分析。首先進(jìn)行焊接過(guò)程的熱傳導(dǎo)計(jì)算，把得到模型的瞬態(tài)溫度數(shù)據(jù)信息存儲(chǔ)在1個(gè)文件中，然后讀入每一時(shí)刻的溫度作為模型的載荷，進(jìn)行焊接模型的應(yīng)力應(yīng)變計(jì)算。由于高斯熱源模型能很好地反映大多數(shù)焊接過(guò)程中的熱流分布特征，其在焊接熱傳導(dǎo)的數(shù)值仿真中應(yīng)用甚廣。高斯熱流分布示意如圖4所示，函數(shù)表達(dá)式為:

圖3 焊縫中線處的圓筒橫剖面Fig.3Cross section of pipe model

式中:q（r）為半徑r處的表面熱流，W/m2;q（O）為熱源中心處的熱流量最大值;c為熱源集中系數(shù);r為距熱源中心的距離。

圖4 高斯熱源熱流分布示意圖Fig.4Heat flux distribution of Gauss heat source

熱傳導(dǎo)計(jì)算中考慮了相變潛熱，取為270 kJ/kg，固相線溫度為1 465℃，液相線溫度為1 544℃;環(huán)境溫度設(shè)為25℃，冷卻條件為空冷，表面換熱系數(shù)取常數(shù)10 W/（m2·℃）;考慮輻射散熱，發(fā)射率取為0.5。加熱過(guò)程中時(shí)間步長(zhǎng)最大取0.5 s，冷卻過(guò)程取自動(dòng)時(shí)間步長(zhǎng)。為了兼顧計(jì)算的效率和精度，在焊縫附近區(qū)域網(wǎng)格劃分得很細(xì)，在遠(yuǎn)離焊縫區(qū)劃分得較粗。沿焊縫方向根據(jù)熱源的位置，在網(wǎng)格重劃分模型中分4步實(shí)現(xiàn)焊接過(guò)程模擬，平板焊接模型的有限元網(wǎng)格如圖5所示。完整網(wǎng)格模型的有限元網(wǎng)格與網(wǎng)格重劃分模型第4步的網(wǎng)格相同，沿縱向均勻劃分，沿橫向由密變疏。圓筒焊接模型的網(wǎng)格劃分方法與重劃分步驟與平板焊接模型相似，在此不贅述。為了簡(jiǎn)化，以FM（Full Model）表示完整模型，RM（Rezoning Model）表示網(wǎng)格重劃分模型。平板完整模型、平板網(wǎng)格重劃分模型、圓筒完整模型、圓筒重劃分模型分別表示為FMⅠ，RMⅠ，F(xiàn)MⅡ，RMⅡ。平板和圓筒焊接過(guò)程中的瞬時(shí)溫度場(chǎng)分布云圖如圖5和圖6所示。

2.2 應(yīng)力應(yīng)變計(jì)算

應(yīng)力應(yīng)變計(jì)算所使用的有限元模型與熱傳導(dǎo)計(jì)算模型的節(jié)點(diǎn)和單元完全相同，單元類型由熱單元改為力學(xué)單元，材料特性由熱物性參數(shù)變?yōu)闊崃W(xué)性能參數(shù)，同時(shí)增加了力學(xué)邊界。在平板有限元模型中只限制了結(jié)構(gòu)的剛體位移，在圓筒有限元模型中設(shè)置了一端簡(jiǎn)支的約束條件。將熱傳導(dǎo)計(jì)算的節(jié)點(diǎn)瞬時(shí)溫度讀入有限元模型，開始應(yīng)力應(yīng)變的增量步計(jì)算。圖7為圓筒焊接后的Mises應(yīng)力分布云圖。

圖7 圓筒焊后Mises應(yīng)力云圖Fig.7Mises Stress distribution of pipe model after welding

3 計(jì)算結(jié)果

3.1 溫度場(chǎng)計(jì)算結(jié)果

平板焊縫長(zhǎng)度中心處橫截面上3點(diǎn)的瞬時(shí)溫度結(jié)果見圖8，測(cè)點(diǎn)分布于平板表面，與焊縫的距離分別為0 mm，20 mm和44 mm。從圖中不難看出，焊接過(guò)程中各點(diǎn)經(jīng)歷了加熱和冷卻2個(gè)階段。距離焊縫越近的點(diǎn)，溫度上升過(guò)程越快，達(dá)到的最高溫度也越高。在冷卻階段，各點(diǎn)的溫度逐漸趨于一致。網(wǎng)格重劃分模型與完整模型的計(jì)算結(jié)果完全一致，說(shuō)明采用網(wǎng)格重劃分技術(shù)在計(jì)算熱傳導(dǎo)過(guò)程時(shí)滿足精度要求。

3.2 變形場(chǎng)計(jì)算結(jié)果

平板焊接模型的變形結(jié)果見圖9，曲線表示了橫向收縮、角變形沿焊縫方向的分布情況。從圖中可以看出，橫向收縮在焊縫長(zhǎng)度中間段有最大值，角變形隨焊縫長(zhǎng)度增加略有變大，但整體上看，2種變形沿焊縫方向都較為均勻。

=圓筒焊接模型中A和B兩點(diǎn)（標(biāo)示于圖3）的瞬態(tài)徑向位移示于圖10。從圖中可以明顯看出，當(dāng)點(diǎn)所在區(qū)域處于加熱階段時(shí)，焊縫金屬發(fā)生膨脹，徑向位移變大，在熱源遠(yuǎn)離該處，焊縫因冷卻而收縮，這時(shí)徑向位移減小。由于A，B兩點(diǎn)在環(huán)焊縫中受熱2次，可在位移圖中看到2個(gè)波峰。

由圖9和圖10可知，平板和圓筒的完整模型與網(wǎng)格重劃分模型計(jì)算結(jié)果相一致，說(shuō)明網(wǎng)格重劃分模型預(yù)測(cè)的焊接變形結(jié)果沒有精度損失。

3.3 應(yīng)力場(chǎng)計(jì)算結(jié)果

平板表面堆焊模型的殘余應(yīng)力結(jié)果以等值線的方式繪制見圖11和圖12。由圖11可見，縱向殘余應(yīng)力的分布情況，焊縫區(qū)域存在較大的拉應(yīng)力，遠(yuǎn)離焊縫區(qū)則呈現(xiàn)較低的壓應(yīng)力以使構(gòu)件保持自身的平衡。圖12為平板焊后的橫向殘余應(yīng)力結(jié)果，在焊縫兩端存在較大的壓應(yīng)力，焊縫區(qū)的橫向應(yīng)力水平較低。圖11和圖12表明，完整計(jì)算模型與網(wǎng)格重劃分模型計(jì)算的焊接殘余應(yīng)力結(jié)果非常吻合。

圓筒焊接模型的0°和180°縱向剖面的殘余應(yīng)力結(jié)果見圖13和14，圖中反映了圓筒焊后內(nèi)外側(cè)的應(yīng)力分布情況。焊縫起始位置處（0°）的環(huán)向應(yīng)力峰值較低，而焊縫中間段（180°）的環(huán)向應(yīng)力峰值則超過(guò)了材料的屈服強(qiáng)度;同時(shí)，圓筒內(nèi)側(cè)的應(yīng)力要高于外側(cè)的應(yīng)力水平。圖14中的結(jié)果表明，圓筒內(nèi)側(cè)的軸向應(yīng)力為拉應(yīng)力，外側(cè)的軸向應(yīng)力為壓應(yīng)力，0°位置處的應(yīng)力峰值低于180°位置處的應(yīng)力水平。從圖中結(jié)果可以看出，完整模型和網(wǎng)格重劃分模型的橫向殘余應(yīng)力和縱向殘余應(yīng)力完全吻合，而且有限元分析預(yù)測(cè)的環(huán)向應(yīng)力和軸向應(yīng)力的分布規(guī)律與Burdekin［6］關(guān)于低碳鋼圓筒環(huán)焊縫的研究結(jié)果相一致。

3.4 有限元計(jì)算時(shí)間

平板與圓筒焊接模型的有限元網(wǎng)格節(jié)點(diǎn)數(shù)和計(jì)算時(shí)間如表3所示，2種模型在使用網(wǎng)格重劃分技術(shù)后，計(jì)算時(shí)間均減少了約30%。由于網(wǎng)格重劃分模型中的節(jié)點(diǎn)數(shù)一般少于完整模型的節(jié)點(diǎn)數(shù)，也即模型中總的自由度數(shù)前者少于后者，內(nèi)存的需求也必將因網(wǎng)格重劃分技術(shù)而得以減少。對(duì)大型結(jié)構(gòu)而言，計(jì)算時(shí)間和內(nèi)存大小是決定數(shù)值模擬能否實(shí)現(xiàn)的主要因素。計(jì)算結(jié)果表明，網(wǎng)格重劃分技術(shù)可以解決這方面的問(wèn)題。本文暫時(shí)沒有研究熱源后端的網(wǎng)格粗化問(wèn)題，若能選擇合理的粗化準(zhǔn)則并優(yōu)化網(wǎng)格設(shè)計(jì)將可能大幅提高計(jì)算效率，同時(shí)保證結(jié)果的可靠性。

4 結(jié)語(yǔ)

1）焊接熱源影響下的區(qū)域具有高度非線性特征，遠(yuǎn)離熱源區(qū)域?yàn)槿醴蔷€性，可以采用網(wǎng)格重劃分技術(shù)進(jìn)行模擬。

2）使焊縫區(qū)細(xì)密網(wǎng)格隨熱源不斷移動(dòng)，與完整的有限元模型相比，減少了節(jié)點(diǎn)自由度，提高了計(jì)算效率并保證了應(yīng)力和變形結(jié)果的精度。

3）通過(guò)網(wǎng)格重劃分模型與完整模型的計(jì)算結(jié)果可知，提出的網(wǎng)格重劃分方法在充分保證結(jié)果精度的前提下，大約減少了30%的計(jì)算時(shí)間。

［1］BROWN SB，SONGH．Rezoninganddynamic substructuring techniques in FEM simulations of welding processes［J］．ASME Journal of Engineering for Industry，1993，155:415－423．

［2］LINDGREN L E，HAGGBLAD H A，et al．Automatic remeshing for three-dimensional finite element simulation 0f welding［J］．Computer Methods in Applied Mechanics and Engineering，1997，147（3）:401－409．

［3］MURTI V，VALLIAPPANS．Numericalinverseisoparametricmappinginremeshingandnodalquantity contouring［J］．Computers and Structures，1986，22（6）: 1011－1021．

［4］HIBBITT D，KARLSONB，SORENSONP．ABAQUS analysisuser'smanual［M］．ABAQUSVersion6.5 Documentation．

［5］WANG Rui，SHERIF R，et al．Numerical and experimental investigations on welding deformation［J］．Trans．JWRI，2008，37（1）:79－90．

［6］BURDEKIN F M．Local stress relief of circumferential butt welds in cylinders［J］．British Welding Journal，1963，10 （9）:483－490．