張玉存 田孟奇 李亞彬
燕山大學電氣工程學院,秦皇島,066004
環(huán)形鍛件是核電、石化等領域的關鍵基礎部件,其成形質量直接影響裝備的總體水平和可靠性,然而鍛件中存在的鍛造缺陷會影響鍛件的成形品質[1]。對含缺陷鍛件的傳熱模型進行研究,實時測量鍛造過程中熱態(tài)鍛件的表面溫度,作為定解條件代入傳熱模型,判斷鍛造過程中內部溫度場的異常變化,及時有效地調整鍛造工序,進而提高鍛件的內部質量和鍛造精度,可為鍛件缺陷的定性識別提供理論支撐,對提高生產達標鍛件的良品率和鍛造的節(jié)能降耗具有重要的現(xiàn)實意義[2]。
近年來,國內外學者對包括環(huán)筒鍛件在內的各種環(huán)形鍛件溫度場傳熱模型進行了大量的研究。范春利等[3]在熱傳導理論的基礎上,針對含內部缺陷的鍛件建立了二維傳熱模型,提出了一種求解得到鍛件表面溫度場以識別試件內部缺陷形狀和位置的計算方法。沈功田等[4]應用紅外熱成像檢測技術,對不同壁厚的碳鋼管道內部和不同深度的壁厚減薄缺陷進行了加溫或降溫過程的紅外熱成像檢測試驗,指出材料熱導率和厚度、溫度激勵方式是影響檢測靈敏度的關鍵因素。關明等[5]、JANIK 等[6]通過試驗模擬實際生產中的鍛件,建立了Mn18Cr18N護環(huán)鋼鍛造中的溫度場,對在一定保溫措施下鍛件的溫度場變化情況進行了描述。崔曉龍等[7]在分析大型鍛件傳熱的基礎上建立了大型鍛件的物理模型,并對其溫度場的變化規(guī)律進行了對比分析,證明了熱處理過程中數(shù)值模擬的可行性。王海亮等[8]針對內部含矩形缺陷試件,建立了二維物理和數(shù)學傳熱模型,通過有限體積法求解了瞬態(tài)溫度場分布。CARLONE等[9]利用有限元仿真了鋼的淬火過程,并成功揭示了溫度場與固-固相變過程的內在關系。沈立華等[10]利用紅外熱像儀獲得外壁溫度,在簡化一維模型下反推內壁溫度,考慮了缺陷內空氣對流換熱及表面輻射對反推內壁溫度的影響。曹春梅等[11]結合導熱反問題,將紅外檢測得到的外壁溫度作為強加邊界條件,反推內壁熱流密度和對流換熱系數(shù)分布的普適解析級數(shù)解。上述研究是在一定理想假設條件下進行的,換熱條件和物理場選擇一維或二維條件下的穩(wěn)態(tài)溫度場,解析模型的缺陷邊界大多為規(guī)則面。在鍛造過程中,鍛件表面的溫度場與內部的結構、熱物性及表面與外界環(huán)境的熱交換是相互影響的。袁龍蔚[12-13]從缺陷體流變學理論上證明,鍛造過程中由于裂縫塑性流變變形會形成溫度場,裂縫處傳熱是連接鍛件內部熱傳導和表面溫度場之間過渡的一部分,是造成溫度場非均勻性的原因之一。
本文通過理論分析方法構建鍛件溫度場傳熱模型,利用相應的數(shù)學物理方法計算溫度場分布函數(shù),通過傳熱模型直接明確地分析和表示裂縫條件下的鍛件溫度場,為鍛造過程中裂縫的識別提供了理論計算方法;有限元仿真應用Comsol代入理論中的推導公式進行數(shù)值計算,驗證理論計算的必要性;試驗中,熱像儀首次測量的表面溫度為初始條件溫度,作為附加定解條件代入傳熱模型,獲得環(huán)形鍛件溫度場任一時刻、任一位置溫度值的解析解,通過對比測量溫度的相對誤差分析試驗精度,試驗將測量溫度近似作為真實溫度;最后,應用MATLAB對裂縫表面溫度值的解析解和測量溫度進行擬合,觀察各自的溫度變化趨勢,進一步驗證傳熱模型的可行性。
導熱微元體始終遵循能量守恒定律,任一時間段均存在熱平衡關系[14-15]。導入微元體總熱流量與微元體內熱源生成熱的和等于導入微元體總熱流量與微元體熱力學能增量的和。如圖1所示,在無內熱源條件下,有
式中,Q為單位時間內傳遞的熱量,即微元體熱力學能增量;ρ、c分別為鍛件的密度和質量熱容;r為微元體半徑。
對微元體的能量收支作熱平衡分析,有
其中,Vv為微元體體積,其質量為 ρVv,dτ內溫度升高所需的熱量為 Φ(r+θdr)drdτ,θ∈(0,1),質量熱容c表示單位質量微元體升高單位溫度所需熱量。
圖1 柱坐標下微元體的導熱、對流熱平衡分析Fig.1 Analysis of thermal and convective heat balance of micro bodies under column coordinates
瞬態(tài)導熱通過微元體的熱流量隨時間改變,式(1)等號右邊第2項為微元體內部存在的熱量蓄積。聯(lián)系傅里葉導熱定律、Robin邊界條件,單位時間內由導熱、對流通過微元曲面dS所傳遞的熱量即熱流量
聯(lián)立各熱力學方程,可得關于微元體的熱平衡關系:
其中,T、Ta分別表示材料溫度和環(huán)境溫度。假如微元體被加熱,即材料溫度高于裂縫處的空氣溫度,則視為熱阻性缺陷。dT dτ為正值,T-Ta為負值,式(2)等號右邊負號是為了對這一事實進行確認。由能量守恒得到的熱平衡關系可知,鍛件與環(huán)境的對流散熱量約等于微元體內能的減少量,也可以通過集中參數(shù)法得出與式(2)基本一致的導熱微分方程式。
應用推導的熱平衡方程對鍛件、裂縫的傳熱進行分析,建立鍛件與裂縫的傳熱微分方程:
定解條件為
式中,Sw*、Sw、SL分別為環(huán)形鍛件內外表面和裂縫邊界面;hw*、hw、hL分別為環(huán)形鍛件內外表面?zhèn)鳠嵯禂?shù)和裂縫表面?zhèn)鳠嵯禂?shù);ρL、cL、VL分別為裂縫高溫干空氣的密度、質量熱容和體積;ρ、c、V分別為鍛件的密度、質量熱容和體積;T、TL、Tf*、Tf分別為材料、裂縫以及環(huán)形鍛件內外表面環(huán)境溫度。
將得到的TL作為瞬態(tài)溫度場的傳熱微分方程的Robin邊界條件,結合實時測得的鍛件表面溫度,建立瞬態(tài)溫度場的傳熱微分方程:
其中,鍛件傳熱微分方程的外表溫度Tw(r2,φi,zi,τ)由紅外熱像儀循環(huán)掃描獲取?;谒矐B(tài)傳熱原理,裂縫對溫度場的影響以傳熱函數(shù)的形式作為鍛件溫度場的Robin條件,考慮熱流連續(xù)性,環(huán)形鍛件裂縫Robin條件如下:
式中,h?r為修正傳熱系數(shù);k?r為修正熱導率裂縫邊界面上溫度的法向方向導數(shù);nx、ny、nz分別為裂縫邊界面上溫度的法向方向余弦;kL為高溫干空氣熱導率,其值見文獻[15];TL分別對 r、φ、z求偏導,是由坐標變換后復合函數(shù)求偏導的鏈式法則生成的。
對于式(6)中的修正傳熱系數(shù),由于導熱的瞬態(tài)熱流量隨時間而變化,故對相鄰單位時間τ~τ+1內微元體的熱流量進行分析:
為保證 l、dφ 有非零解,且 l≠0、dφ≠0 ,從而有
取
則單位面積的熱流量為熱流通量,瞬態(tài)導熱時內部存在熱量蓄積,熱流通量與內部溫度分布隨時間而變化。
于是引入修正熱物性參數(shù)的傳熱微分方程邊界條件:
為簡化參量和計算,取 F=hS(ρcV),引入修正熱物性參數(shù):
代入定解條件式(4),進一步求解裂縫的溫度場函數(shù):
從對流換熱微分方程組出發(fā),聯(lián)系自然對流換熱的準則方程:Nu=C(GrPr)m(L δ)n,計算對流傳熱系數(shù)h=Nuk δ,其中,Nu為努塞爾數(shù),Pr為普朗特數(shù),k為流體熱導率,可得對流傳熱系數(shù)的表達式:
傳熱系數(shù)是與鍛件溫度相關的函數(shù),Th、Td分別表示對流面的高溫函數(shù)和低溫函數(shù),對于環(huán)形鍛件外表面,對應表示外表面溫度Tw、外側環(huán)境溫度Tf;Tm為外表面與外側環(huán)境的平均溫度。
將熱流量作為內外表面溫度的中間變量,分離變量并積分,代入傅里葉導熱定律,可得鍛件溫度與內表面溫度關系:
當半徑r→r2,溫度T→Tw時,可得鍛件內外表面溫度關系:Tw*=Tw-r2(ln r2-ln r1)dT dr。根據上式計算鍛件內側溫度作為Th代入式(8)重新定義熱物性參數(shù)并計算,得到內表面環(huán)境溫度函數(shù):
其中,C、m、n參照努塞爾數(shù) Nu;Pr、v、k的取值見文獻[15]。內側環(huán)境溫度函數(shù)Tw*代入裂縫的溫度場函數(shù)TL,進而使函數(shù)的熱物性參數(shù)變量減小。
最后,為了分離變量以求解偏微分方程,定解條件中的Robin非齊次邊界條件需轉化成齊次形式,考慮其邊值函數(shù)與時間τ無關,作滿足空間變量 r的輔助函數(shù)[16-17]:Tb(r)=(ar+b)Tf+(cr+d)Tf*。輔助函數(shù)的選取有一定任意性,形式、構造方法不唯一,在滿足確定的齊次邊界條件和齊次泛定方程條件下,經計算結果是一致的,代入邊界條件式(7),由對應系數(shù)相等解待定系數(shù)a、b、c、d。令ζ=T-Tb,代入泛定方程并確定相應參數(shù),使得關于ζ的定解問題具有齊次邊界條件和齊次的泛定方程。代入邊界條件式(7),得到關于ζ的齊次邊界條件式及最終的輔助函數(shù):
應用分離變量法求解溫度場傳熱微分方程,各待定系數(shù)、積分常數(shù)可以通過代入裂縫邊界條件、傳熱微分方程定解條件及特征函數(shù)的正交性來確定,迭加分離方程基本解可以得到鍛件溫度場函數(shù)的完全解[18]:
其 中 ,ω2=β2+η2為 假 設 常 量 符 合 關 系 ,N=1 K。分離方程 R(βir)可改寫為v階貝塞爾微分方程的 標 準 解 形式:R(βir)=C1Jv(βir)+C2Yv(βir),Jv(βir)、Yv(βir)分別為第一類、第二類貝塞爾通解函數(shù)。
聯(lián)系初始條件 ζ(r,φ,z)|τ=0=T0-Tb,代入附加邊界條件 T(φi,zi,τ)|r=r2=Tw(τ),可確定式(10)中的系數(shù) Iij,即
其中,Jv+1(βir)、Jv-1(βir)、Yv+1(βir)、Yv-1(βir)是由貝塞爾函數(shù)的微分性質[19]推導的導數(shù)與階次v間的微分關系和遞推公式。聯(lián)系v階貝塞爾微分方程的通解函數(shù)可得對應的通解表達式。代入裂縫邊界條件:
聯(lián)系式(12)可確定 βi值以及最終的含裂縫環(huán)形鍛件溫度場傳熱模型:
應用Comsol有限元軟件建立上文理論設定條件下的計算模型,對熱態(tài)環(huán)形鍛件的瞬態(tài)傳熱進行仿真,獲取含裂縫環(huán)形鍛件的溫度場分布,提取全方向角的軸向多點表面溫度數(shù)據。應用Comsol分析含裂縫環(huán)形鍛件的傳熱性質時,考慮有限元計算需足夠精細化網格,同時保證合理的計算時間,選取如下尺寸的環(huán)形鍛件進行仿真。軸向高度為0.5 m,內外徑分別為0.30 m、0.24 m,裂縫尺寸為深度0.015 m(近表)~0.03 m(環(huán)壁中心)、寬度1~3 mm,高度0.1 m,其他參數(shù)同上文。全方向角、不同軸向高度仿真溫度分布數(shù)據見表1。
表1 環(huán)形鍛件表面的不同高度的圓周上仿真溫度分布數(shù)據Tab.1 The temperature distribution data on the circumference of the different heights of the ring forging surface
由表1可以看出,環(huán)形鍛件上出現(xiàn)了一些溫度低于正常溫度的數(shù)據(用黑體表示),通過選取一維繪圖組對不同尺寸裂縫的溫度數(shù)據進行具體分析。同一時刻條件下,內外環(huán)境溫度不同條件下的外表溫度曲線見圖2??紤]環(huán)形鍛件內外環(huán)境溫差,對于較小或較大尺寸的裂縫,在紅外檢測最小分辨率內的溫度變化幅度更加明顯。
同一時刻下,考慮環(huán)形鍛件內外環(huán)境溫差,對比裂縫邊界不同條件下的外表溫度曲線見圖3??梢钥闯?,環(huán)形鍛件的裂縫邊界為對流時,對于較小或較大尺寸的裂縫,在紅外檢測最小分辨率內的溫度變化幅度更加明顯。
綜上所述,考慮熱態(tài)環(huán)形鍛件內外環(huán)境溫差和裂縫邊界的有限空間自然對流,對于較小或較大尺寸裂縫的識別是必要且可行的。對比圖2、圖3可以得出,毫米數(shù)量級的裂縫寬度具備一定的識別度,裂縫寬度為1 mm時區(qū)分已不明顯,即常規(guī)計算條件下微米數(shù)量級寬度的裂縫已較難識別。為了明確顯示表面溫度變化,應用MATLAB對仿真所得的環(huán)件表面溫度數(shù)據進行簡單的消噪、精簡處理,用環(huán)件表面展開的曲面溫度分布圖表示,見圖4。
圖2 內外環(huán)境溫度不同條件下的外表溫度曲線Fig.2 Internal and external ambient temperature under different conditions of the surface temperature curve
圖3 對比裂縫邊界不同條件下的外表溫度曲線Fig.3 Curve of external surface temperature under different conditions of fracture boundary
圖4 裂縫對應表面的曲面溫度分布圖Fig.4 The surface temperature distribution of thecorresponding surface of the crack
圖4 中,參照xy坐標以及右側反色表,分析缺陷的方位角、軸向高度可知裂縫的大致形態(tài)及位置。鍛件在軸向0.2~0.3 m間溫度出現(xiàn)較大的區(qū)域性波動,且弧度范圍內任意軸向截面大致呈溝狀降低。溫度越高,高溫干空氣熱導率越大,鍛件熱導率越小,高溫干空氣的熱導率遠小于熱態(tài)鍛件的熱導率,為熱阻性缺陷。
為了驗證瞬態(tài)傳熱理論所得瞬態(tài)溫度場傳熱模型的可行性,選取含裂縫熱態(tài)環(huán)形鍛件作為試件,鍛件和裂縫尺寸參數(shù)同上文。試驗中,紅外熱像儀鏡頭與被測鍛件的直線距離為5 m,垂直距離為0.5 m。應用SRJX-8-13箱式電阻爐加熱環(huán)形鍛件,加熱至1 000℃后取出置于滾動臺上,調試紅外熱像儀對鍛件表面溫度進行檢測,鍛件表面的紅外熱圖像見圖5。
圖5 含裂縫熱態(tài)環(huán)形鍛件熱像圖Fig.5 Thermal image with cracks for hot ring forgings
此時,可以準確獲取任一時刻下環(huán)形鍛件表面溫度場,某一時刻的表面3D溫度場見圖6。運用IRBIS 3 professional導出裂縫對應表面的軸向溫度值,選取不同軸向高度、π方向角上解析解與測量溫度、仿真溫度進行比較,見表2。
圖6 鍛件表面3D溫度場分布Fig.6 Forging surface 3D temperature field distribution
表2 不同軸向高度、π方向角上解析解與測量溫度、仿真溫度的比較Tab.2 Different axial height,the direction of the angle analysis and measurement of temperature,simulation temperature comparison
分析表2中數(shù)據可知,解析解與測量溫度、仿真溫度與測量溫度的最大相對誤差分別為0.24%、0.31%,瞬態(tài)傳熱理論計算的傳熱模型解析解滿足精度要求(解析解是將熱像儀首次掃描環(huán)形鍛件表面溫度作為初始溫度信息和其他定解條件代入傳熱模型求解的),驗證了傳熱模型的正確性;基于理論推導公式與邊界條件設定下的Com?sol仿真結果滿足精度要求,驗證了仿真分析的可行性。為了進一步驗證傳熱模型的可行性,應用MATLAB對表2中環(huán)形鍛件表面溫度數(shù)據進行擬合,對比曲線見圖7。圖7中,裂縫邊緣存在局部溫度升高的區(qū)域性波動,符合缺陷體流變學理論。
綜上所述,傳熱模型解析解與試驗測量溫度的變化趨勢大致相同;同時,解析解與試驗測量溫度、仿真溫度與試驗測量溫度的相對誤差滿足精度要求,仿真的有限元計算建立在理論推導公式的基礎之上,說明傳熱模型可以通過熱像儀獲取外表面溫度為初始條件下,實現(xiàn)對熱態(tài)環(huán)形鍛件內裂縫的識別。
圖7 裂縫對應表面的軸向溫度變化Fig.7 Cracks correspond to the axial temperature changes of the surface
相對于一般邊界條件設定常量的處理方法,考慮熱態(tài)環(huán)形鍛件內外環(huán)境溫差和裂縫邊界的有限空間自然對流,對于較小或較大尺寸的裂縫都具有更好的識別度,證明了理論計算的必要性。試驗中,傳熱模型解析解、仿真溫度分別與試驗數(shù)據的對比結果滿足相對誤差的精度要求,驗證了傳熱模型和仿真分析的可行性。
[1] 任運來,聶紹珉,牛龍江,等.大型鍛件內部空洞缺陷修復條件[J].機械工程學報,2008,44(2):248-252.REN Yunlai,NIE Shaomin,NIU Longjiang,et al.Res?toration Conditions of Hollow Voids in Large Forgings[J].Journal of Mechanical Engineering,2008,44(2):248-252.
[2] 戴文遠,李長有.紅外定量檢測缺陷的熱傳導理論分析模型[J].紅外,2014,34(4):43-46.DAI Wenyuan,LI Changyou.We Theoretical Model of Heat Conduction Theory for Infrared Quantitative De?tection Defects[J].Infrared,2014,34(4):43-46.
[3] 范春利,孫豐瑞,楊立.基于紅外測溫試件內部缺陷的識別算法研究[J].工程熱物理學報,2007,28(2):304-306.FAN Chunli,SUN Fengrui,YANG Li.Study on Recog?nition Algorithm of Internal Defects Based on Infrared Thermometry[J].Journal of Engineering Thermophys?ics,2007,28(2):304-306.
[4] 沈功田,李濤,姚澤華.高溫壓力管道紅外熱成像檢測技術[J].無損檢測,2002,24(11):473-477.SHEN Gongtian,LI Tao,YAO Zehua.Infrared Ther?mal Imaging Detection Technology for High Tempera?ture Pressure Pipes[J].Nondestructive Testing,2002,24(11):473-477.
[5] 關明,付赟秋,常志梁,等.大鍛件鍛造過程中溫度場測定及其結果分析[J].鍛壓技術,2012,37(2):6-9.GUAN Ming,F(xiàn)U Yunqiu,CHANG Zhiliang,et al.De?termination of Temperature Field and Its Results in Forging Process of Large Forgings[J].Forging Technol?ogy,2012,37(2):6-9.
[6] JANIK M,DYJA H.Modelling of Three-dimensional Temperature Field Inside the Mould during Continu?ous Casting of Steel[J].International Journal of Heat and Mass Transfer,2004,20:177-182.
[7] 崔曉龍,萬妮麗.大型鍛件熱處理過程的數(shù)值模擬研究[J].熱處理,2005,20(4):12-16.CUI Xiaolong.WAN Nili.Study on Numerical Simula?tion of Heat Treatment Process for Large Forgings[J].Heat Treatment,2005,20(4):12-16.
[8] 王海亮,范春利,孫豐瑞,等.二維內部缺陷的紅外瞬態(tài)定量識別算法[J].紅外與激光工程,2012,41(7):1714-1720.WANG Hailiang,F(xiàn)AN Chunli,SUN Fengrui,et al.In?frared Transient Quantitative Recognition Algorithm for Two-dimensional Internal Defects[J].Infrared and Laser Engineering,2012,41(7):1714-1720.
[9] CARLONE P,PALAZZO G S,PASQUINO R.Finite Element Analysis of the Steel Quenching Process:Temperature Field and Solid-solid Phase Change[J].Computers&Mathematics with Applications,2010,59(1):585-594.
[10] 沈立華,范春利,楊立,等.缺陷對紅外測溫反推內壁溫度的影響分析[J].紅外技術,2005,27(3):250-253.SHEN Lihua,F(xiàn)AN Chunli,YANG Li,et al.Influence of Defects on Inverted Temperature Measurement by Inferred Wall Temperature[J].Infrared Technology,2005,27(3):250-253.
[11] 曹春梅.圓筒內壁熱流密度和對流換熱系數(shù)的紅外熱診斷研究[J].激光與紅外,2007,37(9):849-851.CAO Chunmei.Infrared Thermal Diagnosis of Heat Flux and Convection Heat Transfer Coefficient in Cyl?inder Wall[J].Laser and Infrared,2007,37(9):849-851.
[12] 袁龍蔚.帶缺陷流變性材料裂尖斷裂過程區(qū)的熱力學性和電磁性[J].湘潭大學自然科學學報,1997,19(3):29-37.YUAN Longyu.Thermodynamic and Electromagnetic Properties of Fracture Tip Fracture Zone in a Rheolog?ical Material with Defects[J].Natural Science Jour?nal of Xiangtan University,1997,19(3):29-37.
[13] 袁龍蔚.缺陷體流變學——Ⅱ.裂縫萌生區(qū)的平衡與演化以及缺陷體位移場[J].湘潭大學自然科學學報,1997,19(3):39-55.YUAN Longyu.Rheology of Defective Bodies—II.Balance and Evolution of Crack Initiation Zone and Displacement Field of Defect[J].Natural Science Journal of Xiangtan University,1997,19(3):39-55.
[14] 涂虬.熱平衡法在傳熱學教學中的應用[J].武漢工程職業(yè)技術學院學報,2000,12(1):63-67.TU Qiu.Application of Thermal Balance Method in Teaching of Heat Transfer[J].Journal of Wuhan En?gineering Institute,2000,12(1):63-67.
[15] 陶文銓.傳熱學[M].西安:西北工業(yè)大學出版社,2006:272-277.TAO Wenquan.Heat Transfer[M].Xi’an:Northwest?ern Polytechnical University Press,2006:272-277.
[16] 嚴鎮(zhèn)軍.數(shù)學物理方程[M].合肥:中國科學技術大學出版社,1996:68-74.YAN Zhenjun.Equations of Mathematical Physics[M].Hefei:University of Science&Technology Chi?na Press,1996:68-74.
[17] 傅立葉.熱的解析理論[M].桂質亮,譯.北京:北京大學出版社,2008.FU Liye.Light Translation.Analytical Theory of Heat[M].Gui Zhiliang,Trans.Beijing:Peking University Press,2008.
[18] 張玉存,張雷強,付獻斌.熱態(tài)軸類含空洞鍛件內部溫度場分布[J].中國機械工程,2015,26(19):2672-2676.ZHANG Yucun,ZHANG Leiqiang,F(xiàn)U Xianbin.Analyses of Interior Temperature Field Distribution for Hot Axial Forgings with Void[J].China Mechani?cal Engineering,2015,26(19):2672-2676.
[19] 克羅克斯頓.數(shù)學物理方程導論[M].載安英,譯.北京:高等教育出版社,1982:275-282.Croxton.Introduction to Equations of Mathematical Physics[M].Zai Anying,Trans.Beijing:Higher Edu?cation Press,1982:275-282.