考慮固相變的鈦合金焊接殘余應力數(shù)值仿真方法

酈 羽，萬正權，唐文勇

（1.中國船舶科學研究中心，江蘇無錫 214082；2.上海交通大學船舶海洋與建筑工程學院，上海 200011）

0 引 言

焊接溫度場、殘余應力和變形的實驗測量既耗時又昂貴，而且還可能存在與測量相關的不確定性[1]。因此，焊接計算模擬經(jīng)常被用來預測不同類型材料和接頭在焊接過程中和焊接后的瞬態(tài)殘余應力場和變形的分布和大小。由于焊接過程中的溫度梯度非常大，焊縫在小體積內(nèi)是由復雜的非均勻組織組成。在冷卻過程中焊縫可能會發(fā)生固態(tài)相變，相關的體積變化和相的具體材料特性會對焊接殘余應力和變形產(chǎn)生顯著影響。許多焊接模型忽略了相變的影響，F(xiàn)erro等[2]認為這種影響是不可忽略的。

少數(shù)研究人員研究了Ti-6Al-4V 的相變動力學，其中大多數(shù)研究是針對Ti-6Al-4V 焊接的具體的案例，忽略了復雜相變過程影響，因為Ti-6Al-4V 隨化學成分、溫度歷史和相變前的初始相形態(tài)而變化。Malinov 等[3]利用電阻率技術研究了等溫條件下的相變動力學，并利用差示掃描熱法研究了不同冷卻速率下、連續(xù)冷卻條件下的相變動力學，利用JMAK 方程[4]和疊加性規(guī)則確定了β相轉化為α相的分數(shù)；Ahmed 和Rack[5]研究了在β相轉變點以上的冷卻速度對冷卻過程中相變的影響，并確定在高于410 ℃/s 的高冷卻速度下形成馬氏體；而Gil Mur 研究了回火熱處理過程中馬氏體分解為α和β的過程，并在更低的冷卻速率下觀察到完全馬氏體組織；Kelly Kampe 等[6]模擬了Ti-6Al-4V 在激光金屬沉積過程中的微觀結構演變；Fan 使用Malinov 等[3]的JMAK 等溫轉變動力學參數(shù)，數(shù)值研究了相變對Ti-6Al-4V 激光成形的影響；Elmer 等[7]利用同步x 射線衍射技術實驗測量了Ti-6Al-4V 氣體鎢極弧焊（GTAW）過程中的轉化動力學；Longuet 等[8]開發(fā)了Ti-6Al-4V的通用多相模型，并應用于直接激光制造和激光焊接工藝。

本文分別采用JMAK 方程和K-M 方程表達鈦合金固態(tài)相變過程中的擴散轉變和非擴散轉變建立Ti-6Al-4V 焊接過程中固態(tài)相變的數(shù)學模型，確定數(shù)學模型中的參數(shù)求解方法，并基于ABAQUS 軟件開發(fā)考慮固態(tài)相變的熱-冶金-力學耦合的有限元計算方法，該方法可以考慮溫度與組織變化對焊接殘余應力的影響。通過一個對接接頭的激光焊接的有限元算例，驗證所建立的固態(tài)相變的Ti-6Al-4V 焊接有限元計算模型的可行性。該工作對于定量分析Ti-6Al-4V 焊接殘余應力的影響因素具有參考意義。

1 焊接殘余應力數(shù)值計算方法

焊接是一個涉及許多復雜物理現(xiàn)象的過程，焊接過程中焊件的冶金、熱和力學變化之間存在多場相互作用，包括：焊接熱循環(huán)引起的熱應力以及伴隨的熔化、凝固或固態(tài)相變，局部膨脹引起的轉換應力，非彈性變形產(chǎn)生的熱，由應力或應變加速的相變、相變潛熱。對于解決工程實際問題的計算方法而言，為了兼顧計算效率與計算精度，在有限元模型中往往重點考慮主要因素，有時要忽略次要因素。

1.1 溫度場的計算

由于焊接過程的物理復雜性，采用一種簡化方法僅基于熱方程的求解，其中熱源是模擬的，而不是規(guī)定焊接模型的體積熱流輸入。由傅里葉定律和能量守恒推導出的熱方程如式（1）所示，用于求解瞬態(tài)溫度場（T）在時間（t）和空間（x，y，z）上的溫度場：、，

式中，T是溫度，ρ為密度，Cp為比熱，k為導熱系數(shù)，Q?V為體積熱通量（W/m3）。為了求解熱方程，必須指定熱導率、密度和比熱容。對于經(jīng)歷相變的材料，相變潛熱和微觀結構演變也必須考慮。

采用三維錐形高斯熱源模型對激光焊接過程進行熱分析。圖1 所示為三維錐形模型，其徑向功率密度分布為高斯分布，軸向為線性分布，可以得到更精確的深熔激光束焊接結果[9]。圓錐熱源的上表面熱流最大，下表面熱流最小，如式（2）所示。

圖1 對接焊縫采用的三維高斯型熱源模型Fig.1 Three dimensional Gaussian heat source model for butt weld

式中，r0為特定高度z處熱源的半徑，q0為最大體積功率密度，r為當前內(nèi)部點的半徑，

分布參數(shù)r0從圓錐區(qū)域的頂部到底部線性減小，如式（4）所示。

式中，ze是錐體區(qū)域的上表面，zi是下表面。通過對物體上的體積熱流Q進行積分，并將其重新作為總熱量輸入，可以得到最大的體積功率密度，如式（5）所示。

式中，H=ze-zi，h=z-zi，e 為自然對數(shù)的底數(shù)，η為熱源效率。熱源效率或能量傳遞效率定義為被輻照材料吸收的能量與激光輸出功率的比值。

焊接過程中，工件與周圍空氣之間通過對流和輻射的方式進行熱交換。式（6）和式（7）分別將表面對流和輻射造成的熱損失定義為邊界條件。

式中，T為當前溫度，T0為環(huán)境溫度，Tabs為絕對零度溫度，ε是發(fā)射率或發(fā)射熱輻射的能力，而σ為Stefan-Boltzmann 常數(shù)（5.68×10-8J/K4?m2?s）。輻射的影響在較高的溫度下（靠近焊縫）更大，而在較低的溫度下（遠離焊縫）不顯著。Yang等[10]指出：Ti-6Al-4V 的發(fā)射率在760 ℃以下的溫度下幾乎不受溫度的影響，其發(fā)射率在0.3 以下；而在760 ℃以上，它開始氧化，其發(fā)射率隨著溫度的升高而增加；在1000 ℃以上，其發(fā)射率幾乎達到1.0。相比之下，對流效應在溫度較低時更顯著，而在溫度較高時則不顯著。用于傳熱系數(shù)和輻射常數(shù)的校準值為hconv=10，ε=0.4。

1.2 組織計算

鈦合金在中高溫區(qū)的相變可以被看成是兩個均析轉變（一個發(fā)生在加熱過程，另一個發(fā)生在冷卻過程），通過過往研究，鈦合金α→β和β→α相變被認為是兩個互逆的擴散型相變，擴散轉變可以用如公式（8）所示的Johnson-Mehl-Avrami方程（簡稱JMAK方程）[4]來描述。

式中，f為時刻t時產(chǎn)物相的體積分數(shù)，feq是產(chǎn)物相T溫度下的平衡體積分數(shù)，與溫度相關的k( )T和與溫度無關的n是定義相轉變動力學的JMAK方程參數(shù)。

為了將等溫模型的應用擴展到焊接過程中的非等溫相變，使用了Scheil加法規(guī)則[11]。Scheil規(guī)則將瞬時小的等溫時間增量步長疊加的總和近似表示任意連續(xù)的溫度變化。通過設置一個較小的計算時間步，就可以采用階梯狀的等溫變化曲線來近似焊接過程的連續(xù)變化曲線。

對于每一個計算時間步，引入虛擬時間t*0，虛擬時間的表達式為

式中，t*0是在溫度T1的等溫轉變過程中達到f(t0,T0)所需的虛擬時間，feq(T1)是產(chǎn)物相在溫度T1下的相平衡分數(shù)。

在快速冷卻過程中，β相轉變?yōu)轳R氏體相的過程為非擴散型轉變，采用Koisten-Marburger 方程（簡稱K-M方程）進行描述。

式中，fM是馬氏體相的體積分數(shù)，kM是依賴材料的參數(shù)，MS是馬氏體開始溫度。

1.3 應力場計算

計算焊接應力場時，將溫度場的計算結果以熱載荷的方式加載到應力計算模型中，并將組織（相）的計算結果代入材料的參數(shù)模型中，采用熱-彈塑性有限元方法進行力學分析。假設材料的彈塑性行為具有各向同性硬化規(guī)律，鈦合金焊接過程中的總應變包括以下幾部分：

式中，εTotal是總應變，εe是彈性應變，εp是塑性應變，εth是熱應變，εc是蠕變應變，εvp是粘塑性應變，εtp是相變塑性應變?？倯冇缮鲜鰬兎至拷M成。假定蠕變、粘塑性和突變分量引起的應變可忽略不計，因此不包括在計算中。由于忽略了粘塑性效應，假定屈服應力與應變率無關，且取決于塑性應變和溫度。根據(jù)Lindgren 的研究[12]，焊接模擬必須考慮彈性、塑性和熱應變來計算殘余應力。熱應變對焊接誘導殘余應力的影響最大，因此，熱膨脹系數(shù)是應力計算中最重要的材料參數(shù)之一。Lindgren[13]進一步指出，焊接冷卻階段是導致殘余應力和變形的重要階段。在材料經(jīng)歷相變的情況下，所產(chǎn)生的體積變化可以包括在熱應變分量中。

2 Ti-6Al-4V固態(tài)相變模型的建立

2.1 加熱過程相變數(shù)學模型建立

在加熱階段，Ti-6Al-4V 經(jīng)歷了從α相到β相的快速擴散控制轉變，其中包括β的形核和生長以及α/β界面的運動導致的α的溶解，這是由β穩(wěn)定劑在界面上的傳輸引起的。β的體積分數(shù)在β轉變溫度（Tβ）以上達到100%，并在熔化溫度下保持穩(wěn)定。假設β生長增加量是通過相等量的α溶解，該過程為擴散轉變，采用JMAK 方程來描述。為了避免球狀、魏氏組織、網(wǎng)籃組織和晶界α相等相形態(tài)的差異，假設了一個簡化的單α相，如公式（13）所示，計算加熱過程中β相的體積分數(shù)。

式中，t*0是在溫度T1的等溫轉變過程中達到fβ(t0,T0)所需的虛擬時間，feqβ是β相的平衡分數(shù)。

閱讀不是讓學生記住多少知識，而是在閱讀的過程中最終學會閱讀。因此，家長在閱讀繪本的過程中可以滲透相關繪本閱讀的策略，比如，看圖獲取信息并進行推斷，對后來發(fā)生的事進行預測，閱讀中的自我提問，結合學生的實際等。家長們可以根據(jù)自己的優(yōu)勢和特長自主設計大膽嘗試，讓孩子們充分體驗閱讀帶給他們的快樂，進一步提升他們的閱讀能力，培養(yǎng)他們熱愛閱讀的濃厚興趣和良好習慣，從而讓孩子們逐漸學會閱讀繪本，構建閱讀力。

2.2 冷卻過程相變數(shù)學模型建立

在冷卻階段，冷卻速率決定了轉變動力學和形成的相。在低于β相變點的溫度下，當冷卻速率小于20 K/s時，β通過成核和擴散反應分解成α[5]。采用JMAK方程描述。

式中，t*0是在溫度T1的等溫轉變過程中達到fα(t0,T0)所需的虛擬時間，feq

α是α相的平衡分數(shù)。

在β相變點以上溫度快速冷卻，將β轉變成另外一種形式的α相，即細針狀馬氏體相α'，其不同于平衡α相，由于相鄰片層和界面的結晶取向差異較大，且重復率較低，因此具有較高殘余應力[14]。

根據(jù)Ahmed和Rack的研究[5]，對于大于410 ℃/s的快速冷卻速率，將形成完全馬氏體組織，并抑制向α相的擴散轉變。

對于20 ℃/s到410 ℃/s之間的中等冷卻速率，從β相變點以上溫度開始冷卻，將導致β向晶界的塊狀α'部分轉變，和靠近原β晶界的晶體內(nèi)馬氏體αm的部分轉變。將α'納入αm中，考慮到缺乏有關該問題的精確數(shù)據(jù)，可進行另一種選擇，并從可用于轉變的β相的量中減去β相的當前平衡分數(shù)feqβ，以防止完全轉變?yōu)轳R氏體α'。

式中，fαm是αm的分數(shù)，fβeq是β的平衡分數(shù)。在這種情況下，從可用于轉化的β相中減去β的當前平衡分數(shù)，以防止完全轉化為α'。

2.3 相變數(shù)學模型中的參數(shù)確定

為了求解JMAK 方程的產(chǎn)物相，需要確定方程中的參數(shù)k( )T和n，以及產(chǎn)物相的平衡相，圖2 顯示了JMatPro 模型確定的α相和β相的平衡分數(shù)。由于實際的母材組織中存在β相，在低溫條件下，β相體積分數(shù)預測不能為0。

圖2 α相和β相平衡相分數(shù)隨溫度變化的模型Fig.2 Equilibrium fraction of α and β phases varying with temperature

本文基于Szkliniarz 和Smotka[15]的試驗數(shù)據(jù)，在不同加熱和冷卻速度下材料處于不同溫度時的α、β兩相體積分數(shù)，對其兩組試驗數(shù)據(jù)，對方程（8）進行改寫，采用最小二乘法擬合求解。

通過擬合的數(shù)據(jù)發(fā)現(xiàn)，n0基本不隨溫度變化，k0（Ti）隨溫度而變化。為了使JMAK 方程計算的值更準確，將擬合的k0（Ti）和n0作為初值，將該溫度下的JMAK 方程計算值fJMAK(t,Ti)與該溫度下的試驗測得值fmeasure(t,Ti)之差的平方（fJMAK（t,Ti）-fmeasure（t,Ti））2，與試驗測得值fmeasure(t,Ti)的比值記為誤差值Eerr，通過搜索求解k(Ti)和n的值，使得Eerr值最小而作為該溫度下的一組最優(yōu)解（k（Ti）,n）。

為了較為準確地表達出k（T）隨溫度的變化，Elmer 等將隨溫度變化的參數(shù)k（T）表示為式（21）所示的形式，p和Q均為不隨溫度變化的系數(shù)，R為氣體常數(shù)。將以上優(yōu)化的一組不同溫度下的k（T）值，通過公式（21）擬合求得p、Q值。

對于K-M 方程中的參數(shù)求解，通過GLEEBLE 熱模擬試驗機對試樣進行熱膨脹試驗得到熱膨脹曲線。在較大冷卻速度范圍內(nèi)冷卻組織轉變?yōu)轳R氏體相，通過雙切線法確定馬氏體相變的開始溫度MS。由K-M方程的逆運算得到馬氏體相變參數(shù)kM的值，計算公式如下：

圖3 顯示了β的平衡相分數(shù)以及α'的平衡相分數(shù)，它們使用K-M 方程的kM或γ和MS函數(shù)確定。選擇γ等于0.015，而MS等于650°C。

圖3 用Koistenen-Marburger模型計算β和α'的相分數(shù)，作為Ms和γ的函數(shù)Fig.3 Calculation of phase fractions of β and α'with Koistenen-Marburger model as the functions of Ms and γ

3 Ti-6Al-4V的結果與討論

3.1 焊接過程有限元模型的建立

為了驗證所建立的鈦合金焊接過程中的固態(tài)相變過程，本文采用一個50 mm×50 mm×1.5 mm的板進行對接焊計算，激光焊接速度為50 mm/s。計算模型的網(wǎng)格尺寸為0.5×0.5×0.5 mm3，單元數(shù)量為60 000。在熱分析中，采用2 階三維20 節(jié)點二次擴散傳熱六面體單元（DC3D20）。在力學分析中，采用了一階三維8結點線性六面體簡化積分單元（C3D8R）。力學分析中位移的有限元形狀函數(shù)通常要比熱分析高一個維度。簡化積分單元，C3D8R 使用低階積分來形成單元剛度，并提供更好的近似真實行為。Ti-6Al-4V隨溫度變化的熱物理材料性能參數(shù)如圖4所示。

圖4 Ti-6Al-4V隨溫度變化的熱物理材料特性Fig.4 Thermo-physical material properties of Ti-6Al-4V varying with temperature

3.2 焊接過程中相的成分變化

由第2 章所建立的Ti6Al4V 相變模型可知，Ti6Al4V 的相變轉變類型主要與焊接過程中的三個變量相關：當前溫度、冷卻/加熱速率和峰值溫度。通過對焊接過程中這三個變量的判斷，可以根據(jù)相變方程計算出相應的相產(chǎn)物體積分數(shù)。為了實現(xiàn)焊接過程的相變模擬，在熱分析過程中，本文使用ABAQUS 的USDFLD 子程序的內(nèi)部狀態(tài)變量作為冷卻/加熱速率和峰值溫度的函數(shù)，計算不同微觀結構相的相變體積分數(shù)。在力學分析中，通過UEXPAN 子程序調(diào)用USDFLD 子程序的相體積分數(shù)場變量，將不同相的熱膨脹特性參數(shù)與相變體積分數(shù)疊加運算，從而考慮焊接相變對熱應變影響，最終得到殘余應力的影響。

在模型計算過程中，當前溫度為場變量NT11（℃），冷卻/加熱速率定義為狀態(tài)變量SDV3（℃/s），峰值溫度定義為狀態(tài)變量SDV4（℃），場變量FV1 定義為β相，場變量FV2 定義為α相，場變量FV3 定義為馬氏體相αm。圖5～7 為焊接過程中，0.25 s、1.0 s、4.0 s 三個時刻仿真計算場變量變化；圖8～10 為焊接過程中，0.25 s、1.0 s、4.0 s 三個時刻鈦合金板的組成成分體積分數(shù)變化。Ti6Al4V 試樣板初始組成成分：β相體積分數(shù)為0.075 26，α相體積分數(shù)為0.924 74，馬氏體αm相體積分數(shù)為0。

圖5 焊接過程中0.25 s時刻仿真計算場變量變化Fig.5 Simulation calculation of field variable change at 0.25 s during welding

圖8 焊接過程中0.25 s時刻組織成分體積分數(shù)變化Fig.8 Volume fraction changes of different micro-phases at 0.25 s during welding

由峰值溫度SDV4 可知，焊接過程中，只在焊縫及熱影響區(qū)域的一個較小范圍內(nèi)峰值溫度超過650 ℃，也即焊接過程中的相變主要發(fā)生在該區(qū)域。由當前溫度NT11 可知，隨著熱源的移動，溫度場由熔池（溫度超過1650 ℃的區(qū)域）向外逐漸降低，呈橢圓分布。溫度場的分布特征決定了冷卻/加熱速率SDV3 的分布，由圖5 和圖6 可知，隨著熱源移動，由內(nèi)到外呈現(xiàn)冷卻速率到加熱速率的過渡，冷卻速率由內(nèi)到外逐漸減小，最內(nèi)部超過-410 ℃/s 的區(qū)域呈長橢圓分布，逐漸過渡到-410～-20 ℃/s 的區(qū)域。

圖6 焊接過程中1.0 s時刻仿真計算場變量變化Fig.6 Simulation calculation of field variable change at 1.0 s during welding

對圖5進行分析，在0.25 s時刻，焊縫及熱影響區(qū)域當前溫度都在980 ℃以上，所以都轉化為β相，與圖8 計算的相變組成相吻合。在1.0 s時刻，整個焊縫及熱影響區(qū)域的溫度在650 ℃以上，所以馬氏體αm相不會產(chǎn)生，冷卻速率在-410～-20 ℃/s 區(qū)域時，生成α相，與圖9 的相變組成相吻合。在4.0 s 時刻，也即冷卻過程，焊縫及熱影響區(qū)域周圍有相當大部分區(qū)域的溫度已經(jīng)降到300～650 ℃之間，焊縫及熱影響區(qū)域的冷卻速率在-410～-20 ℃/s，因此該階段焊縫及熱影響區(qū)會產(chǎn)生馬氏體αm相，由圖10可以發(fā)現(xiàn)此時焊縫及熱影響處主要是馬氏體αm相和α相，含有少量β相，與圖7也吻合得較好。因此，仿真計算過程基本實現(xiàn)了第2章中Ti6Al4V固態(tài)相變的相變模型過程。

圖7 焊接過程中4.0 s時刻仿真計算場變量變化Fig.7 Simulation calculation of field variable change at 4.0 s during welding

圖9 焊接過程中1.0 s時刻組織成分體積分數(shù)變化Fig.9 Volume fraction changes of different micro-phases at 1.0 s during welding

圖10 焊接過程中4.0 s時刻組織成分體積分數(shù)變化Fig.10 Volume fraction changes of different micro-phases at 4.0 s during welding

圖11 顯示了焊接熱循環(huán)過程中，焊縫中心距起始位置20 mm 處積分點位置α相、β相和α'相體積含量變化的數(shù)值模擬結果。初始組織由約0.924 74α和0.075 26β組成。在加熱過程中，隨著溫度的迅速升高，α 相分數(shù)迅速下降到0，而β相增加到1，實際上，在熔融溫度以上也應該變成0，但為了簡單起見，在β轉變點溫度以上，它仍然維持在1。當材料開始冷卻到β相變點溫度以下時，根據(jù)冷卻速率的不同，一個不同的微觀結構由β相的擴散或非擴散轉化產(chǎn)生。由圖11 可知，冷卻過程中β相先轉變?yōu)棣料?，冷卻到一定溫度后，α相的體積分數(shù)保持不變，剩下的β相轉變?yōu)轳R氏體αm相，達到β相在該溫度下的平衡相后，各相成分保持不變。在冷卻至室溫過程中，最終相主要由0.472 7α，0.488 64αm以及0.038 66β相組成。

圖11 焊接中心處α、β和α'相分數(shù)隨溫度變化的時間-溫度歷史Fig.11 α,β and α'phase fractions at the weld center line chang?ing with temperature in the time-temperature history

3.3 焊后殘余應力和平面變形

從圖12 和圖13 可以看出，在數(shù)值焊接模型中引入相變效應會對焊接過程結束時焊縫周圍殘余應力的最終分布產(chǎn)生一定改變。相變模型的y向應力分量S22 峰值比無相變模型大100 MPa，且焊縫寬度上的應力梯度更陡，殘余應力范圍更大。同樣，通過相變計算得到的x向應力S11 和z向應力S33 值也會增大，且改變了拉壓應力的分布情況，但只是增大了很小的一部分。具有相變效應的變形結果比沒有相變的變形結果略有增加，彎曲變形和角變形的幅度在0.1 mm以內(nèi)，由于模型尺寸較小，變形不是很明顯。

圖12 中橫剖面處三個主應力方向的殘余應力分布Fig.12 Residual stress distribution in three principal stress directions at the middle cross section

圖13 工件焊后殘余應力和變形的分布Fig.13 Distribution of residual stress and deformation of workpiece after welding

基于以上仿真結果分析，由于α和β單個晶體的尺寸匹配良好，這意味著在新相的形核和溶解過程中，沒有任何原因改變初始晶向，但是鈦合金材料在不同相之間的體積比存在較大差異。因此，Ti-6Al-4V 中的相變會引起焊縫體積變化，從而導致焊接殘余應力的一定變化，但是不會顯著改變殘余應力的分布情況。

4 結 語

本文建立了Ti-6Al-4V 焊接過程中固態(tài)相變的數(shù)學模型，確定了數(shù)學模型中的參數(shù)求解方法，并基于ABAQUS 軟件開發(fā)了考慮固態(tài)相變的熱-冶金-力學耦合的有限元計算方法?；谠撚嬎隳Ｐ停瑢i-6Al-4V 對接板激光焊接進行了計算驗證，得出了焊接接頭的組織包含α相、馬氏體αm相和少量β相。并基于以上結果，進行了焊接殘余應力計算。通過對計算結果分析可以發(fā)現(xiàn)，考慮相變后殘余應力y向分量S22增加明顯，x向分量S11和z向分量S33的分布規(guī)律有較明顯改變，但是幅值較小。鈦合金材料在不同相之間的體積比存在較大差異，因此相變會引起焊接殘余應力的一定變化，但是不會顯著改變殘余應力的分布情況。