,,,,
(1.上海交通大學上海市激光制造與材料改性重點實驗室,上海200240;2.高新船舶與深海開發裝備協同創新中心,上海 200240;3.創斯達科技集團(中國)有限責任公司,江蘇南通 226300)
焊接是金屬材料加工制造中一種重要的材料連接方式,連接可靠性非常高,加工周期較短、加工效率高、易于實現自動化,在工業制造領域的應用非常廣泛,是現代制造業中不可或缺的材料連接方法。然而焊接過程相當復雜,涉及到金屬材料局部的快速熔化及材料的快速凝固,局部溫度的快速升高及冷卻導致焊接過程中和焊接后工件內部產生相當大的焊接殘余應力以及焊接變形。而焊接過程中產生的殘余應力及變形影響著工件的精度和使用性能,降低工件的使用壽命。
此外,在焊接溫度場和焊接應力的測量方面,由于焊接過程及之后的冷卻凝固過程速度非常快,導致焊接溫度場和應力應變的測量與分析相當困難。焊接數值模擬方法的提出與使用為測量分析焊接殘余應力和變形提供了一定的參考,特別是結構件焊接過程模擬為焊接工藝及參數、焊接順序等的制定給予一定的指導。同時,現代計算機技術的發展推動了焊接數值模擬技術的發展,使得研究工作者的工作越來越深入,新方法不斷涌現,對焊接方面的實際生產和科研起到了推動作用[1]。
有限元法是古典變分法與分片插值法相結合的產物。它先將全域根據實際情況分成很多個小單元,在單元范圍內用低次多項式分片插值,再將它們進行組合,使其能形成全域內的函數,再用以逼近問題的真解。這樣不但能避免古典方法尋找基函數的困難,而且不規則剖分與差分方法相比靈活性和適應性更大,應用范圍極廣,能計算物理和工程中的各種復雜問題。
在焊接中常用的熱彈塑性有限元法集合了傳熱有限元分析過程和結構彈塑性有限元分析過程兩個部分,能夠較好地反映材料的結構參數、焊接工藝參數以及不同的焊接順序等因素對焊接殘余應力和焊接變形等方面的影響[2]。
在Ronsenthal-Rykalin解析式[3]中,熱源模型可根據焊接結構件的厚度和尺寸形狀將熱源簡化為點、線、面三種形狀的熱流簡化模型。對于大型焊接結構件,焊接時熱量的傳導沿著3個方向,因此可以將熱源看成點狀熱源,關系式為

式中 Q為熱源中心溫度,Q=ηUI;η為焊接熱源的效率系數;α為熱擴散率;D為距焊接熱源中心的距離。
而對于細棒的對接,溫度場在細棒的截面上均勻分布,可以將細棒的焊接截面上的各點看成是均溫的,因此可以將該焊接熱源看成面熱源,關系式為

該計算方法以集中熱源為基礎,假定熱物理性能參數不變,不考慮相變與結晶潛熱對焊接殘余應力應變的影響,將焊接工件的幾何形狀簡單地歸為無限長、大、薄,這種模型對于遠離焊縫區的計算結果較為精確,但是對于熱影響區誤差較大。該模型的計算結果雖然精度較差,但是計算過程較為簡單,因此在工程上得到廣泛應用。
為了提高計算結果的精度,研究者們參考焊接時的焊接熱源分布形狀,又提出了高斯熱源模型、半球狀熱源模型、橢球狀熱源模型及雙橢球狀熱源模型來對應不同情況下的焊接模型。
焊接時的焊接熱源是通過一定的作用面積將熱量傳遞到焊接工件其他部位,在這個作用面積上熱量分布并不均勻,分布規律為中心部位熱量最高,邊緣部位熱量較低。費利德曼將這個作用面積上的熱量分布近似地用高斯數學模型進行描述,如圖1所示。數學表達式為


圖1 高斯熱源的數學模型
式中 r為模型上某一點距離熱源中心處的距離;q為距離焊接熱源中心r處的熱源密度;Qm為熱源中心處的最大熱源密度;R為焊接熱源的有效加熱半徑。焊條電弧焊、鎢極氬弧焊等焊接電弧挺度較小、溫度梯度相對較小的焊接方法,運用高斯熱源模型能夠取得較精確的計算結果。
與高斯熱源模型相比,半球狀熱源模型和橢球狀熱源模型更能適應電弧挺度、溫度梯度相對較大的焊接方法,例如高能束焊接。另外考慮到運用激光焊等方法焊接時,焊接熔池并不是球對稱的,對于這種情況橢球狀熱源模型較半球狀熱源模型更為合適。
半球狀熱源模型的數學表達式為

橢球狀熱源模型數學表達式為

式中 a,b,c分別為在x,y,z軸方向上的半軸長。
加拿大學者J.Goldak[4]提出的雙橢球熱源模型是根據橢球熱源模型改進而來。在橢球熱源模型中,橢球前半部分的溫度梯度變化沒有實際焊接中劇烈,并且橢球后半部分的溫度梯度變化較為平緩,不滿足溫度梯度和電弧挺度較大的激光焊等焊接特點,因此提出了雙橢球熱源模型。在雙橢球熱源模型中,將半側的橢球分為2個1/4的橢球,前半部分是一個1/4的橢球,后半部分是一個1/4的橢球。前半部分分配的能量系數為n1,后半部分分配的能量系數為n2,n1和n2滿足的條件為:n1+n2=1。
前半部分橢球熱源模型的數學關系式為

后半部分橢球熱源模型的數學關系式為

式中 a,b,c可取不同的值。
通常的焊接方法采用高斯熱源;對于電弧沖力較大的熔化極氬弧焊、激光焊接采用雙橢球的效果較好,還可以將熱源分為兩部分,采用高斯分布熱源作為表面熱源,雙橢球分布作為內熱源[5]。
除了上述常用熱源模型外,研究者們又根據模型的實際情況提出了基于標準幾何形狀的其他熱源模型,使其能夠在特定的焊接工藝中取得較好的計算結果。根據計算結果發現,熱源的選定對整個模型的計算過程及計算結果有著至關重要的作用,需根據實際情況選取合適的熱源模型[6]。
在焊接模擬中,針對瞬態移動熱源計算時間過長、計算效率過低的問題,清華大學的鹿安理和蔡志鵬等人[7]提出了高斯分段熱源。在高斯分段熱源模型中,將焊縫分成若干段,每一段對應一個時間步,每一段內各個點同時施加高斯熱量密度,段狀熱源模型如圖2所示。

圖2 高斯段狀熱源數學簡易模型
這樣既減少了點狀熱源在焊縫上移動所需時間步,大大縮短了計算時間,提高了計算效率,又通過合適的分段數量,保證了模型的計算精度。單位時間高斯熱源在作用面上輸入熱量為

在如圖3所示的焊接長度下,高斯段狀熱源的總輸入熱量Q的數學表達式為


圖3 段狀高斯熱源數學解析式示意
清華大學的王煜等人[8]將該方法等效地應用于電子束焊的雙橢球狀熱源模型中,計算焊接的殘余應力應變,既保證了計算精度,又減少了計算時間,提高了計算效率。
上海交通大學陳震等人[9]認為蔡志鵬等人的分段熱源模型中輸入的熱能密度不隨時間變化,總的熱能輸入量需要對瞬態移動熱源的總輸入量進行計算求得,因此對其進行改進,特點為:①模型中的各段熱能密度輸入隨著時間變化;②分段熱源模型中的熱能輸入過程與瞬態熱源模型一致;③分段熱源模型與瞬態熱源模型之間可以進行簡單地切換,不用再通過復雜的數學積分變換。
自適應網格方法是指計算中在某些變化較為劇烈的區域,如大變形、接觸間斷面和滑移面等,網格在迭代過程中不斷調節,細化網格,做到網格點分布與計算解的調和,從而提高解的精度和效率的一種技術。自適應網格希望在物理解變動較大區域網格自動密集,而在物理解變化平緩區域網格相對稀疏,這樣在保持高效率計算的同時可得到高精度的解。自適應網格技術主要有移動網格方法和局部細化或粗化網格方法。近三十年來,自適應網格方法成為網格方法研究的熱點問題,在一些領域得到廣泛應用。
法國的P.Duranton等人[10]在316L不銹鋼的焊接中運用SYSWELD軟件進行計算,網格劃分方法為自適應網格法,并根據熱源模型的特點和熱源模型的移動軌跡對網格進行“運動中的重新劃分”,如圖4所示。

圖4 自適應網格法網格劃分示意
自適應網格法的提出使得很多復雜工藝的模擬得到實現。在P.Duranton的實驗中,自適應網格法得到的結果與實驗結果一致,并且計算時間是完全細化網格法所消耗時間的20%。
瑞典的Lindgren L E等人[11]在電子束焊中分別運用h-adaptive自適應網格法和非自適應網格法進行對比計算分析,結果發現使用兩種方法得到的計算結果基本相同,但是h-adaptive自適應網格法所需要的計算時間比非自適應網格法縮短了大約60%。
克羅地亞的Mato Peric等人[12]計算不銹鋼EN 10025-2:S355JR的T型接頭MAG焊時分別運用Shell/3D Solid混合單元法和實體單元法進行建模,建模方式如圖5和圖6所示。

圖5 實體單元網格劃分示意

圖6 Shell/3D Solid混合單元網格劃分示意
Mato Peric等人在殼單元和三維實體單元的交界處用邊界約束方程進行約束,約束方程為

式中 Nx、Ny、Nz為法線向邊界的方向余弦值;hc和hr分別為對流換熱系數和輻射換熱系數;qs為邊界換熱系數;Tr為輻射溫度;T∞為室溫。
同時,對三維實體單元和混合單元的計算結果進行對比,結果發現,溫度場的分布在Shell/3DSolid單元交界處溫度梯度過渡平緩,較好地反映出焊接溫度場分布,且殘余應力應變分布也與純實體單元模型得到的結果相似。并且采用Shell/3DSolid混合單元法縮短了42%的計算時間。
將一個大型的復雜結構劃分為若干個子結構,先分別確定各子結構的剛度特性,然后將子結構裝配成整體結構,最后確定整體結構的剛度特性,這種結構分析方法稱為子結構法。采用子結構分析可將一個大型問題化為若干個小問題,將大型的聯立方程組分解為若干組小型的方程組,從而減小計算機的內存,實現微機解大題的作用。
某些結構中,材料的形狀及應力的變化只集中在局部,其他部位幾乎不變,因此將不變部位的結構劃分為若干個子結構,變化部位劃分為其他的子結構。當結構變化時,只需重新形成變化部分的子結構的剛度矩陣,而不必重新形成全結構的剛度矩陣,從而提高計算效率。
Y.G.Duan等人[13]采用有限元軟件SYSWELD計算大型結構件時運用了熱彈塑性有限元子結構法,認為焊接引起的塑性應變和金相組織應定位在焊縫附近,可以采用局部三維模型進行確定。對于整體結構,認為整體變形是由于局部的變形所導致,局部塑性應變發生后,作為初始應變然后轉移到整體模型中,最后得到整體結構的計算結果。
日本學者提出和應用了固有應變ε*的概念,認為焊接固有應變是物體在焊接后存在的塑性應變、溫度應變和相變應變三者之和。也可以表征物體在無外力和內應力的情況下,物體從焊接狀態恢復到自由狀態后,總應變ε減去彈性應變εe的值,即ε*=ε-εe。材料的焊接固有應變是一個隨著材料屬性、材料尺寸及焊接參數變化的量,在使用焊接固有應變值時需先對它進行測量分析[14]。
日本九州工業大學的Toshio Terasaki等人[15]選用不同的焊接參數對不同厚度的低碳鋼板進行平板堆焊的焊接實驗,通過測量焊后試樣的角變形及橫向收縮的大小,得出角變形及橫向收縮與焊接線能量和低碳鋼板厚度的關系曲線,如圖7、圖8所示。
上海交通大學的魏良武博士[16]通過焊接變形實驗和熱彈塑性有限元法測量材料的固有應變,定量分析固有應變與主要影響因素的關系。在SM400A鋼中,得到關于橫向固有應變系數ξ、縱向固有應變K與輸入熱量的關系式。
參數為:板寬300 mm,板長300 mm,板厚1.5~40 mm,熱量200~500 J/mm,焊接速度10 mm/s。
其關系式為

圖7 熱輸入參數對角變形的影響

圖8 熱輸入參數對橫向收縮的影響

重慶交通大學的張繼翔教授等人[17]分析T型焊接接頭的雙側同步焊接固有應變分布,得到在不同的參數下固有應變的變化趨勢為:
(1)在焊縫方向上的橫向收縮的變化趨勢是先增大后減小,且在焊縫的兩端由于焊縫端部效應,產生的橫向固有應變高于焊縫其他部位。
(2)在焊縫方向上的縱向固有應變的變化趨勢為沿焊縫方向由兩端向中間逐漸減小。
(3)T型接頭的腹板厚度和翼板厚度對固有應變的影響較大,這是由于隨著板厚的增加,材料對焊接變形的抗性增大,因此變形減小,固有應變減小。
(4)隨著焊接速度的增大,固有應變隨之減小,這是由于隨著焊接速度的增大,焊接線能量減小導致焊接變形減小,因此固有應變減小。
(5)腹板寬度和翼板寬度對材料焊接固有應變的影響較小。
重慶交通大學的張繼翔教授運用固有應變有限元法計算分析大型橋梁鋼箱梁焊接,驗證了固有應變有限元法的可行性。
合肥工業大學的李萌盛[18]運用初應變法和初應力法研究低碳鋼T型焊接接頭,發現固有應變有限元法中的初應力法和初應變法所得結果與熱彈塑性有限元法所得結果基本吻合,運用固有應變有限元法可以計算大型薄壁結構件的焊接變形,并且大大減少計算時間,提高計算效率,而且在運用固有應變有限元法時可以運用殼單元代替實體單元進行計算,進一步縮短計算時間。
相較以前的經驗公式等方法,計算機數值模擬方法能提供更有力的理論指導,計算機的發展與普及也為數值模擬提供了更方便的物質條件。國內外的研究者們針對焊接計算模擬方面的難點及不足之處提出并運用了很多優化方法,除了分段熱源法、子結構法、固有應變法和混合單元法,還有平行計算法等其他常用方法,使得數值模擬技術為工業生產與科學研究提供更有力的指導。
但是,現在的焊接過程模擬工作還存在很多問題需要解決,比如:(1)高溫下材料的熱物理性能參數的缺乏。這是由于實驗設備方面的限制,無法徹底解決,并且材料熱力學參數的確定過程也較為復雜,得到的參數可能會與材料本身的參數存在一定誤差。(2)焊接電弧的有效功率系數難以確定。(3)材料焊接時的相變過程難以處理。
[1]謝元峰,肖漢斌.基于ANSYS的焊接溫度場和應力的數值模擬研究[D].武漢:武漢理工大學,2006.
[2]商躍進.有限元原理與ANSYS應用指南[M].北京:清華大學出版社有限公司,2005.
[3]Teng T L,Fung C P,Chang P H,et al.Analysis of residual stresses and distortions in T-joint fillet welds[J].International Journal of Pressure Vessels and Piping,2001,78(8):523-538.
[4]Goldak J,Chakravarti A,Bibby M.A new finite element model for welding heat sources[J].Metallurgical and Materials Transactions B,1984,15(2):299-305.
[5]陳家權,肖順湖,楊新彥,等.焊接過程數值模擬熱源模型的研究進展[J].裝備制造技術,2005(3):10-14.
[6]高耀東,張福寬,高俊萍.利用不同焊接熱源模型對結構件進行焊接模擬[J].內蒙古科技大學學報,2013(1):59-63.
[7]蔡志鵬,趙海燕,鹿安理,等.焊接數值模擬中分段移動熱源模型的建立及應用[J].中國機械工程,2002,13(3):208-210.
[8]吳甦,趙海燕,王煜,等.高能束焊接數值模擬中的新型熱源模型[J].焊接學報,2004,25(1):91-94.
[9]沈濟超.大型船體結構焊接變形熱彈塑性有限元數值模擬方法研究[D].上海:上海交通大學,2015.
[10]Duranton P,Devaux J,Robin V,et al.3D modelling of multipass welding of a 316L stainless steel pipe[J].Journal of Materials Processing Technology,2004(153):457-463.
[11]Lindgren L E,Haggblad H A,McDill J M J,et al.Automatic remeshing for three-dimensional finite element simulation of welding[J].Computer Methods in Applied Mechanics and Engineering,1997,147(3-4):401-409.
[12]Peric M,Tonkovic Z,Rodic A,et al.Numerical analysis and experimental investigation of welding residual stresses and distortions in a T-joint fillet weld[J].Materials&Design,2014(53):1052-1063.
[13]Duan Y G,Vincent Y,Boitout F,et al.Prediction of welding residual distortions of large structures using a local/global approach[J].Journal of mechanical science and technology,2007,21(10):1700-1706.
[14]羅宇,鄧德安,江曉玲,等.熱變形的固有應變預測法及實例[J].焊接學報,2006,27(5):17-20.
[15]Terasaki T.Effect of Welding Conditions on Residual Stress Distributions and Welding Deformations in Welded Structure Materials[C]//The Twelfth International Offshore and Polar Engineering Conference.International Society of Offshore and Polar Engineers,2002.
[16]魏良武.固有應變法預測焊接變形的研究及其工程應用[D].上海:上海交通大學,2004.
[17]張繼祥,劉紫陽,王智祥,等.橋梁鋼結構T型接頭雙側同步焊有限元建模[J].電焊機,2014,44(8):77-83.
[18]陳重鈞,李萌盛,王成.固有應變法在T型接頭焊接變形研究中的應用[J].現代焊接,2012(12):16-18.