專利名稱：一種預測熱軋過程板帶溫度場的有限元方法
技術領域：
本發(fā)明屬于軋制技術領域，特別涉及一種預測熱軋過程板帶溫度場的有限元方法。
背景技術：
目前在板帶熱軋過程中，板帶溫度預測方法主要有現場測量、數學模型、有限差分法和有限元軟件仿真。現場測量雖然便捷，但是受主觀因素較大，并且大多只能獲得板帶表面溫度，不能得到全面準確的溫度信息；數學模型極大的依賴于現場試驗條件，計算方法不夠靈活，計算精度往往較低；有限差分法是一種直接將微分問題變?yōu)榇鷶祮栴}的近似數值解法，其理論上的缺陷使得計算精度的提高受到了一定限制；通用有限元軟件由于需要占用大量的內存，降低了計算效率，并且由于其專業(yè)性不強，在求解軋制過程溫度時精度不高。這些因素一定程度上影響了熱軋過程中板帶溫度的準確預測，進而不能夠準確設定軋制過程參數，影響軋制過程的順利進行。
發(fā)明內容
針對現有技術中預測熱軋過程板帶溫度場的各種方法精度低和專業(yè)性不強等缺點，本發(fā)明提供一種熱軋過程預測板帶溫度場有限元法，其目的是提高板帶溫度場預測精度，提高計算效率。
實現本發(fā)明目的技術解決方案如下 利用有限元法求解熱軋板帶溫度場時進行如下假設 軋向尺寸遠遠大于寬度和厚度方向尺寸，因此忽略軋制方向熱傳導 寬度方向上邊界換熱條件和幾何形狀對稱，考慮二分之一斷面 ①采集軋制過程數據，包括軋制參數，材料熱物性參數，單元劃分信息，初始信息，控制參數 軋制參數道次軋制的出口厚度，軋制壓力，軋輥溫度，除磷水流量 材料熱物性參數熱傳導系數，比熱，密度 單元劃分信息寬度單元數和厚度單元數 初始信息軋件溫度，軋件寬度，軋件厚度 控制參數迭代次數 ②對橫斷面進行單元劃分，建立有限元分析模型如圖1所示，對單元節(jié)點進行編號，計算節(jié)點坐標。
根據單元劃分數據、軋件寬度和厚度尺寸建立有限元分析模型，單元和節(jié)點編號沿厚度方向和寬度方向逐漸增加，圖1中，i為單元編號，j為節(jié)點編號，H為厚度，W為寬度。AB和AD絕熱邊界，BC和CD換熱邊界。以A點坐標為零，計算各節(jié)點坐標，在寬度方向和厚度方向上單元均勻劃分。
③根據軋制過程實際條件和軋制階段確定換熱系數h和內熱源強度
熱軋過程包括空冷階段，除磷階段，軋制階段。整個軋制過程不同階段的換熱系數和內熱源強度計算如下 (1)熱軋板帶在空冷過程中，
值為零；其表面換熱方式主要為輻射和自然對流，換熱系數h通過式(1)和式(2)計算 HR＝σ·ε·(T+Tair)(T2+Tair2)(1) 其中HR為輻射系數，σ為Stefan-Boltzman常數，σ＝5.67×10-8W/(m2·K4)；ε為黑度系數，ε與溫度的關系式為ε＝0.125(T/1000)2-0.38(T/1000)+1.1。
熱軋板帶在空冷過程的對流方式為自然對流，其表達式為
其中T(K)為板帶表面溫度；Tair為環(huán)境溫度；b表示板帶寬度。
(2)熱軋板帶在高壓水除磷過程中，
值為零；表面換熱方式為強迫對流和側面輻射，換熱系數h通過式(1)和式(3)計算 HCW＝124.7×w0.663×10-0.00147(T-273.16)(3) 其中w(L/min·m2)為水流密度；T(K)板帶表面溫度。
(3)在軋制過程中，板帶與軋輥之間接觸換熱是主要熱損失方式，忽略塑性變形和摩擦做功，接觸換熱系數與軋制壓力有關。換熱系數h通過式(4)計算 IHTC＝695pm-34400(W/m2K)(4) 式中pm(MPa)-軋制壓力 ④利用有限元基本原理，計算四邊形等參單元的型函數、B矩陣和雅克比矩陣J和雅克比矩陣行列式|J|。
⑤利用空間域有限元離散和時間域有限差分相結合法建立溫度場有限元求解的線性方程組。
(1)以熱力學第一定律為依據建立熱傳導微分方程，假設材料導熱各向同性，二維熱傳導的基本方程為 其中T瞬時溫度(K) ρ材料密度(kg/m3) c材料比熱(J/kg·K)) t時間(s) k熱傳導系數(W/(m·K))
內熱源強度(J/m3) (2)利用歐拉方程在給定邊界條件和初始條件下將二維熱傳導問題變?yōu)榈刃Х汉磉_式求極值問題 每個單元的等效泛函表達式表示為 根據熱傳導問題的變分原理，對泛函式(6)求一階偏導數并置零，根據常用的有限元組合方法，對離散單元進行組裝，把單元的剛度矩陣裝配成整體剛度矩陣，得到二維有限元法求解溫度場的系統(tǒng)方程 其中[KT]-溫度剛度矩陣，[K3]-變溫矩陣，{p}-常數項列式，{T}-溫度列式；E-單元總數；上標e表示每個單元。
對每個單元來說，剛度矩陣、變溫矩陣和常數項可以通過式(8)求解 其中k熱傳導系數(W/(m·K))；ρ材料密度(kg/m3)；c材料比熱(J/(kg·K))；h-換熱系數，N型函數；i，j節(jié)點編號。
(3)利用二點向后差分格式，將系統(tǒng)方程轉化為瞬態(tài)溫度場求解的線性方程組。將系統(tǒng)方程(7)中的溫度對時間偏導數表示為二點向后差分格式 將時間向后差分格式(9)帶入系統(tǒng)方程得到溫度場求解的線性方程組 ⑥對形成的線性方程組(11)進行求解，獲得解為溫度分布。采用一維變帶寬存儲法求解大型線性方程組，從而節(jié)約存儲空間，提高計算效率。一維變帶寬存儲方法可以說明如下 對一般的線性方程組AX＝B，如果系數矩陣A為對稱稀疏矩陣?？梢詮牡谝涣蟹橇阍亻_始進行存儲，直到對角線元素為止
存儲為一維矩陣為 A＝a11 a12 a22 a23 a33 a14 a24 0 a44 a35 0 a55 … an-1，n ann 其中a-為矩陣元素，下表第一項表示行數，第二項表示列數。
⑦根據迭代次數判斷軋制過程中某一軋制階段的溫度計算是否結束。如果此階段溫度計算沒有結束，迭代次數增加，繼續(xù)計算；如果此階段溫度計算結束，那么進行下一步操作，判斷整個軋制過程溫度計算是否結束。
⑧根據軋制過程的總時間判斷整個熱連軋過程溫度計算是否結束。如果整個熱連軋過程溫度計算沒有結束，那么采集數據繼續(xù)進行下一階段計算，如果結束，程序預測溫度過程結束。
有限元求解溫度場的計算程序流程圖如2所示。
本發(fā)明的最大效果是能夠獲得很高的溫度預測精度和整個熱軋過程板帶溫度場分布的詳細信息，為軋制過程提供設定和優(yōu)化參數；本方法應用性強，縮短了計算時間，提高了計算效率。本發(fā)明適用于熱連軋過程的粗軋機和精軋機。
圖1本發(fā)明方法的有限元模型圖， 圖2本發(fā)明方法的有限元法軟件流程圖， 圖3本發(fā)明粗軋過程溫度變化圖， 圖4本發(fā)明精軋過程溫度變化圖， 圖5本發(fā)明開始軋制前溫度分布圖， 圖6本發(fā)明粗軋結束溫度分布圖， 圖7本發(fā)明精軋結束溫度分布圖， 圖1中i為單元編號，j為節(jié)點編號，H為厚度，W為寬度，1為換熱邊界，A為開軋前溫度測量值，B為粗軋結束后溫度測量值，C為平均溫度，D為上表面溫度，E為下表面溫度，F為側面溫度，G為心部溫度。
具體實施例方式 選擇一鋼種熱軋過程為計算對象，利用有限元法分析從出加熱爐到精軋結束整個過程的溫度演變規(guī)律。
某一鋼廠熱連軋過程。粗軋過程有4個機架，第二機架為可逆機架；精軋過程有7個機架?？赡鏅C架軋制5道次，其余機架軋制一道次。
例計算條件見表1。
表1 整個熱連軋過程，共軋制15道次，每一道次軋制的出口厚度見表2所示。
表2 采用表1計算條件和表2道次出口厚度對某一鋼廠熱軋過程板帶溫度演變規(guī)律進行計算分析，計算結果如圖3和圖4，其中3為粗軋過程溫度變化，4為精軋過程溫度變化，圖5，圖6和圖7為不同測量位置的溫度分布等值線。
對此熱連軋過程進行現場溫度測量與跟蹤，開軋表面溫度約為1150℃，粗軋結束溫度值約為990-1000℃，采用有限元法開發(fā)的板帶熱軋溫度預測計算模型，計算的開軋溫度和粗軋結束溫度值與實測值吻合良好，誤差在1％之內，證明了本發(fā)明方法的準確性和可靠性。
本發(fā)明方法針對板帶熱軋過程可以獲得整個軋制過程板帶任何位置的溫度分布，信息詳盡準確。本發(fā)明方法整個熱軋過程計算耗用CPU時間約384ms。本方法預測熱軋板帶溫度極大的提高了溫度預測精度和效率，能夠更好的優(yōu)化熱軋過程參數。
權利要求
1.一種預測熱軋過程板帶溫度場的有限元方法，其特征包括以下步驟
(1)采集軋制過程數據，包括軋制參數，材料熱物性參數，單元劃分信息，初始信息，控制參數
軋制參數道次軋制的出口厚度，軋制壓力，軋輥溫度，除磷水流量
材料熱物性參數熱傳導系數，比熱，密度
單元劃分信息寬度單元數和厚度單元數
初始信息軋件溫度，軋件寬度，軋件厚度
控制參數迭代次數；
(2)對橫斷面進行單元劃分，建立有限元分析模型，對單元節(jié)點進行編號，計算節(jié)點坐標；
(3)根據不同軋制過程，確定邊界換熱系數和內熱源強度
熱軋板帶在空冷過程中，其表面換熱方式為輻射和自然對流，輻射系數表述為
HR＝σ·ε·(T+Tair)(T2+Tair2)
式中HR為輻射系數
σ＝5.67×10-8W/(m2·K4)
ε＝0.125(T/1000)2-0.38(T/1000)+1.1
熱軋板帶在空冷過程的對流方式為自然對流，其表達式為
式中
T(K)為板帶表面溫度
Tair為環(huán)境溫度
b表示板帶寬度
熱軋板帶在高壓水除磷過程中，換熱方式為強迫對流和側面輻射，輻射計算方法同上，對流系數表達式為
HCW＝124.7×w0.663×10-0.00147(T-273.16)
式中
w(L/min·m2)為水流密度
T(K)為板帶表面溫度
在軋制過程中，板帶與軋輥之間接觸換熱是熱損失方式，接觸換熱系數表達式為
IHTC＝695pm-34400(W/m2K)
式中pm(MPa)-軋制壓力；
(4)利用有限元基本原理，計算四邊形等參單元的型函數、B矩陣和雅克比矩陣J和雅克比矩陣行列式|J|；
(5)對有限單元的溫度剛度矩陣和變溫矩陣進行組裝，形成溫度求解的系統(tǒng)方程式
式中
[KT]-溫度剛度矩陣
[K3]-變溫矩陣
{p}-常數項列式
(6)采用一維變帶寬存儲法求解線性方程組，獲得瞬態(tài)溫度場。
全文摘要
一種預測熱軋過程板帶溫度場的有限元方法，屬于軋制技術領域，本方法包括以下步驟(1)采集軋制過程數據；(2)對橫斷面進行單元劃分，建立有限元分析模型，對單元節(jié)點進行編號，計算節(jié)點坐標；(3)根據不同軋制過程，確定邊界換熱系數和內熱源強度(4)利用有限元基本原理，計算四邊形等參單元的型函數、B矩陣和雅克比矩陣J和雅克比矩陣行列式|J|；(5)對有限單元的溫度剛度矩陣和變溫矩陣進行組裝；(6)采用一維變帶寬存儲法求解線性方程組，獲得瞬態(tài)溫度場。本發(fā)明的優(yōu)點能獲得很高的溫度預測精度和整個熱軋過程板帶溫度分布的詳細信息，為軋制過程提供設定和優(yōu)化參數；應用性強，縮短了計算時間，提高了計算效率。
文檔編號G06F17/50GK101178746SQ20071015898
公開日2008年5月14日 申請日期2007年12月18日 優(yōu)先權日2007年12月18日
發(fā)明者李長生, 梅瑞斌, 劉相華 申請人:東北大學