李 斌, 陳 豐, 史良宵
(華北電力大學 能源動力與機械工程學院,河北保定 071003)
?
鍋筒截面瞬態溫度場的導熱正反問題耦合解法
李斌,陳豐,史良宵
(華北電力大學 能源動力與機械工程學院,河北保定 071003)
針對外壁受熱的增壓鍋爐鍋筒,提出了求解其截面瞬態溫度場的導熱正反問題耦合解法.根據鍋筒外壁是否受熱,將其外壁劃分為受熱和不受熱2個區域.對于不受熱區域,沿外壁周向布置熱電偶,根據實際測量所得溫度,利用導熱反問題解法求解該區域的溫度場;對于受熱區域,利用導熱正問題解法求解其溫度場;對于耦合邊界區域,將不受熱區域反問題解法得到的交接邊界處溫度賦值給受熱區域正問題解法作為已知邊界條件,從而實現正反問題耦合,得到整個鍋筒截面的瞬態溫度場.利用Ansys軟件對鍋爐冷態啟動過程中鍋筒的溫度場進行了計算,并與正反問題耦合解法的計算結果進行了對比.結果表明:正反問題耦合解法具有較高的精度,在復雜邊界條件下具有很好的適應性,能夠滿足工程應用的需要.
鍋爐鍋筒; 溫度場; 導熱反問題; 導熱正問題; 正反問題耦合解法
鍋筒是鍋爐最大的厚壁承壓部件,在啟、停過程中及變負荷運行時,承受著由于溫度分布不均勻而產生的熱應力和由內部壓力引起的機械應力的共同作用.頻繁的應力變化會引起鍋筒的疲勞壽命損耗.因此,對鍋筒應力進行在線監測,在鍋爐工作過程中控制鍋筒應力在允許范圍內,進而減小鍋筒的疲勞壽命損耗,這對鍋爐的安全運行具有重要意義[1-5].
溫度場計算是鍋筒應力計算和疲勞壽命分析的基礎,精度將直接影響疲勞壽命分析的結果.傳統的鍋筒溫度場求解根據鍋筒內部換熱條件,在已知結構參數、熱物性參數、初始條件和邊界條件的前提下,通過導熱微分方程得到溫度場.這種方法又稱為導熱問題的直接解法或正問題解法[6-7].該方法根據邊界條件的不同,可適用于復雜邊界條件情況(不同邊界條件分別作用或者幾種邊界條件同時作用).由于鍋筒內部流動和換熱過程復雜,內壁換熱的邊界條件難以準確確定,鍋筒內壁與飽和水和飽合水蒸氣的對流傳熱系數hw和hs通常采用經驗數據,從而導致溫度場計算存在誤差.
為了彌補導熱正問題解法的不足,有學者提出了溫度場計算的導熱反問題解法[8-10].其思路是在鍋筒的外壁選定適當的位置布置熱電偶,然后根據熱電偶測得的溫度值,采用控制容積法,建立控制容積的能量守衡方程,進而求解其瞬態溫度場.Murio[11]根據有限差分法采用反問題解法求解一維導熱反問題;Duda和Taler采用有限差分法將反問題解法運用到形狀簡單的控制體,求解其二維導熱反問題,然后又采用控制容積法和反問題解法求解復雜形狀控制體的多維導熱反問題;Taler等[12]采用反問題解法通過控制容積法的理論對多維導熱反問題進行求解,并通過一維和二維導熱反問題進行精度驗證;Taler等[13]采用反問題解法求解瞬態溫度場并應用于電廠實際在線監測系統.
導熱反問題解法適用于外壁不受熱鍋筒的溫度場求解,可以減小由于內壁傳熱系數采用經驗公式計算帶來的誤差.但對于外壁受熱的情況,如增壓鍋爐鍋筒[14],其上部區域近似絕熱,但其下部外壁部分邊界受爐內火焰輻射換熱,且為了維持微正壓,夾層處熱空氣與鍋筒外壁進行對流換熱,換熱情況復雜,爐內高溫環境下無法布置熱電偶,這種情況不再適宜采用導熱反問題解法.
為此,筆者提出了導熱正反問題耦合解法,將導熱正問題解法與反問題解法相結合,求解復雜邊界條件下某增壓鍋爐鍋筒截面的瞬態溫度場,并利用Ansys軟件對該方法的計算精度進行了驗證.
圖1為某增壓鍋爐的鍋筒,其中鍋筒內徑r1=0.65 m,鍋筒外徑r5=0.74 m.鍋筒內壁換熱條件為上部空間的汽和下部空間的水與鍋筒對流換熱.外壁不同區域的換熱條件不同,圖1中上部白色區域外壁為絕熱條件,下部灰色區域為爐內火焰與鍋筒外壁的直接輻射換熱,兩端陰影部分為夾層處熱空氣與鍋筒外壁的對流換熱.
根據鍋筒外壁的換熱條件,將鍋筒劃分為2部分(見圖2):外壁不受熱區域(采用導熱反問題解法求解)和外壁受熱區域(采用導熱正問題解法求解).

圖1 鍋筒截面簡化模型

圖2 鍋筒截面求解區域劃分示意圖
由于鍋筒上部白色區域外壁絕熱,可以在其邊界上布置熱電偶,通過測得的外壁溫度值,采用導熱反問題解法求得該區域的溫度場.下部灰色區域根據其內外壁邊界條件和初始條件采用正問題解法求解導熱微分方程,得到其溫度場.兩區域交接邊界處采用邊界耦合的方式,將反問題解法求得的兩側耦合邊界的溫度值作為已知的第一類邊界條件帶入正問題解法,實現整個溫度場的求解.
根據正問題解法和反問題解法的需要,整個鍋筒截面劃分的網格如圖3所示,其中S1和S2表示耦合邊界,其余數字表示節點編號.

圖3 鍋筒截面網格示意圖
2.1反問題解法
如圖4所示,在鍋筒外壁節點11~節點17處布置熱電偶.

圖4 反問題解法局部網格示意圖
根據熱電偶測得的外壁溫度值,考慮到外壁存在一定的散熱量,對于鍋筒外壁每個節點所代表的控制容積,采用熱平衡法和傅里葉定律列出各節點的能量守恒方程.外層節點13~節點15的能量守恒方程如下:
(1)
(2)
(3)
式中:c為材料的比熱容,J/(kg·K);ρ為材料密度,kg/m3;λ為材料導熱系數,W/(m·K);Δφ為控制容積角度變化量,rad;Δr為控制容積徑向長度變化量,m;ti為節點i的溫度,℃;qi為節點i的熱流密度,W/m2.
根據式(1)~式(3)可以求解出中間層節點7~節點9的溫度:
(4)
(5)
(6)
其中,α為熱擴散率(α=λ/ρc),m2/s.
同理,根據中間層節點8的能量守恒方程得到內層節點3的溫度:
(7)
由此通過外層節點的溫度,逐次反推,求解得到內層節點的溫度.同理改變不同外層節點位置,可得到整個內層節點的溫度,從而得到鍋筒截面反問題解法求解區域的瞬態溫度場.
2.2正問題解法
根據鍋筒內部換熱條件,在已知結構參數、熱物性參數、初始條件和邊界條件的前提下,采用Simpler算法[2]求解導熱微分方程,從而得到其溫度場.
區域離散采用內節點法,290個節點分別位于188個子區域的中心,此時子區域即為控制容積.圖5為正問題解法的網格劃分示意圖.
導熱微分方程為
(8)
邊界條件如下:

圖5 正問題解法網格劃分示意圖
絕熱邊界
(9)
內壁邊界
(10)
(11)
夾層對流
(12)
輻射換熱
(13)
初始溫度為
(14)
式中:q為熱流密度,W/m2,爐內火焰與鍋筒內壁換熱處近似為第二類邊界條件;hj為夾層外壁處對流傳熱系數,W/(m2·K);tf為夾層外壁流體溫度,℃;τ為時間,s;τ為時間,s;t0為初始溫度,℃;t∞為流體的飽和溫度,℃.
邊界條件的處理如下:當邊界條件為第二、第三類邊界時,采用附加源項法對邊界進行處理,即把由第二、第三類邊界條件所規定的引入或者導出計算區域的熱量作為與邊界相鄰的控制容積的當量源項,有利于統一模式處理邊界條件.
最后采用交替方向線迭代及塊修正的方法求解離散方程,得到整個求解區域的溫度場.
2.3正反問題耦合解法
結合邊界條件的不同,分別采用正問題解法和反問題解法對不同離散區域進行單獨的溫度場求解.
首先根據已知的反問題解法求解區域外壁的溫度值(熱電偶測溫),采用反問題解法求解該區域的溫度場.兩區域交接邊界處(耦合邊界S1和S2)采用正反問題耦合解法,將反問題解法求得的交接邊界處的溫度值作為正問題解法的已知條件(第一類邊界條件),然后通過插值的方式賦值給正問題解法求解區域,并將其作為已知邊界條件,這樣就完成了2個求解區域的溫度傳遞.然后在已知結構參數、熱物性參數、初始條件和其他邊界條件的前提下,求解得到正問題解法求解區域的溫度場.
當求解結束后,再將正問題解法求解得到的邊界溫度重新帶回,進行迭代計算,直到輸入與輸出偏差在允許范圍內,迭代計算結束,至此將溫度場聯系起來耦合求解,得到整個鍋筒截面的瞬態溫度場,計算流程如圖6所示.

圖6 計算流程圖
為了驗證正反問題耦合解法求解鍋筒溫度場的準確性,利用Ansys軟件對該方法進行驗證.
首先,利用Ansys軟件對鍋爐冷態啟動過程中鍋筒的溫度場進行計算.如圖7所示,計算時模型沿鍋筒徑向劃分為5個節點,沿圓周方向劃分為72個節點.鍋筒邊界條件的取值如下:q=50 000 W/m2,hw=1 000 W/(m2·K),hs=2 000 W/(m2·K),hj=142 W/(m2·K),tf=127 ℃,t0=70 ℃.鍋筒內流體的飽和溫度隨時間的變化關系如下:
(15)
然后,將利用Ansys軟件求解得到的啟動過程各時刻對應反問題解法求解區域的鍋筒外壁溫度作為該區域的已知條件,正問題解法求解區域采用的初始條件、邊界條件等與Ansys軟件計算時相同,重新應用本文方法計算鍋筒啟動過程的溫度場.將2種方法的計算結果進行對比,從而驗證正反問題耦合解法的正確性.

圖7 Ansys軟件中劃分的網格示意圖
3.1節點溫度驗證
根據不同邊界條件,對鍋筒中若干個具有代表性的節點(見圖3)溫度進行分析驗證,對比分析正反問題耦合解法和Ansys軟件2種方法的計算結果,溫度隨時間的變化如圖8~圖14所示.
圖8給出了內壁節點5處2種方法計算結果的對比.該節點位于汽空間的鍋筒內壁,即反問題解法求解區域.因為反問題解法需要已知外壁溫度,所以將Ansys軟件計算得到的外壁溫度帶入反問題解法,求解得到內壁溫度,再與Ansys軟件計算所得內壁溫度進行對比分析.由圖8可以看出,2種方法計算結果的吻合度很高.

圖8 節點5處溫度的對比
圖9給出了節點19和節點24處2種方法計算結果的對比.這2個節點位于汽水交接邊界處,由于汽水混合,加之飽合水蒸氣和飽合水與鍋筒內壁的傳熱系數相差很大,換熱復雜,造成該處附近點波動較大,溫差也較大,3 000 s時內外壁溫差很小,溫度趨于均勻.
圖10給出了節點20和節點25處2種方法計算結果的對比.這2個節點分別位于鍋筒內水空間的內外壁處.整個過程中溫度隨時間的變化率比汽空間處小,與圖9相比,圖10中的換熱強度較小,到3 000 s時內外壁溫差較大.

圖9 節點19和節點24處溫度的對比

圖10 節點20和節點25處溫度的對比
節點21和節點26處由于外壁受到爐內輻射換熱的影響,相對常規電站鍋爐的鍋筒而言,外壁溫度高,外壁向內壁換熱.在開始階段,內外壁溫度隨時間的變化率基本一致,當飽和水溫度不變時,該處外壁溫度繼續升高,而內壁溫度則變化較小(見圖11),導致內外壁溫差逐漸增大,產生較大的熱應力.

圖11 節點21和節點26處溫度的對比
圖12給出了節點22和節點27處2種方法計算結果的對比.這2個節點分別位于正問題解法求解區域的夾層對流換熱處的內外壁.由于夾層處外壁流體溫度比鍋筒內水的溫度先高后低,所以造成該處內外壁溫度相對于外層絕熱區域外壁溫度先高后低,同時內外壁溫差最后也比絕熱區域大.

圖12 節點22和節點27處溫度的對比
節點23和節點28分別位于正問題解法求解區域的汽空間內外壁處.相對于圖10中水空間處的節點溫度而言,由于上部飽和水蒸氣與鍋筒內壁的傳熱系數大,所以曲線斜率大,節點的溫度變化更快;同時內外壁溫差比水空間處大;到達穩定狀態時,該區域內外壁溫差較小,溫度分布均勻(見圖13).

圖13 節點23和節點28處溫度的對比
圖14給出了耦合邊界節點1和節點10處溫度的對比.耦合邊界將反問題解法求解得到的邊界溫度作為已知邊界賦值給正問題解法,求解完后再將邊界溫度反代回去,編制程序反復迭代,得到精確解.該處溫度的變化與上部汽空間區域(見圖13)一致.由圖14可以看出,2種方法的計算結果吻合得很好.
由圖8~圖14還可以看出,2種方法在不同區域的求解結果吻合度都很高,說明所提出的正反問題耦合解法具有較高的精度.
3.2鍋筒截面瞬態溫度場的分布
采用正反問題耦合解法求解鍋筒截面的瞬態溫度場分布,不同時間截面瞬態溫度場的分布如圖15所示.

圖14 節點1和節點10處溫度的對比
圖15(a)給出了τ=500 s時鍋筒截面的等溫線.此時下部熱流作用時間較短,鍋筒截面上下部溫差相對較小,整個截面溫度分布較均勻,下部內外壁溫差也較小.
圖15(b)給出了τ=1 000 s時鍋筒截面的等溫線.此時由于下部熱流作用時間較長,影響范圍增大,下部空間內壁溫度高于上部空間,同時內外壁溫差增大,上下部溫差逐漸增大,熱應力也相應增大.

(a) τ=500 s

(b) τ=1 000 s
Fig.15Transient temperature distribution of the drum section at different times
(1)導熱反問題解法思路簡單明了,求解結果更加精確,但僅適用于外壁不受熱區域;導熱正問題解法適用范圍廣,但計算精度不高.
(2)導熱正反問題耦合解法彌補了單一方法求解存在的局限和不足,提高了計算結果的精度,實用性廣泛,可適用于復雜邊界條件下溫度場的求解分析,具有一定的理論意義和較高的工程應用價值.
(3)采用Ansys軟件對導熱正反問題耦合解法進行驗證,對比分析了不同位置處節點溫度隨時間的變化,結果表明導熱正反問題耦合解法具有較高的精度,在復雜邊界條件下具有很好的適應性,能夠滿足工程應用的需要.
[1]李立人, 陳瑋, 盛建國, 等. 鍋爐受壓元件的高溫蠕變-疲勞壽命設計計算方法[J].動力工程,2009,29 (5):409-416.
LI Liren,CHEN Wei,SHENG Jianguo,etal. Creep-fatigue life design and calculation method for boiler pressure elements under elevated temperature[J]. Journal of Power Engineering, 2009, 29(5):409-416.
[2]管德清,莫江春,張學綸,等.電站鍋爐汽包壽命在線監測系統[J].動力工程,2002,22(6):2044-2047.
GUAN Deqing,MO Jiangchun,ZHANG Xuelun,etal. Life on-line monitoring system for boiler drum of power station[J]. Power Engineering, 2002,22(6):2044-2047.
[3]鄭心偉,孫瑜,王曉軍.增壓鍋爐汽包低周疲勞壽命計算方法研究[J].熱能動力工程,2010,25(2):184-189.
ZHENG Xinwei,SUN Yu,WANG Xiaojun. Study of the methods for calculating the low-cycle fatigue life of a supercharged boiler drum [J]. Journal of Engineering for Thermal Energy and Power, 2010,25(2):184-189.
[4]劉彤,徐鋼,龐力平,等.鍋爐爐內承壓部件的蠕變分析及壽命計算[J].動力工程,2004,24(5):631-635.
LIU Tong,XU Gang,PANG Liping,etal. Creep analysis and life calculation of the pressure components inside boiler[J]. Power Engineering, 2004,24(5):631-635.
[5]HU Wensen,LI Bin,CAO Zidong,etal. An inverse method for online stress monitoring and fatigue life analysis of boiler drums[J]. Journal of Chongqing University (English Edition),2009,8(2):89-96.
[6]賈鴻祥,林友新.調峰機組鍋爐汽包溫度場解析[J].西安交通大學學報,1994,28(6):92-98.
JIA Hongxiang,LIN Youxin.The analysis of temperature field in the steam generator drum wall for cyclic loading units[J].Journal of Xi'an Jiaotong University,1994,28(6):92-98.
[7]趙鐵成,沈月芬,朱國楨.電站鍋爐鍋筒溫度場計算——三維非穩態變物性材料不均勻導熱問題有限元分析[J].中國電機工程學報,1997,17(4):217-220.
ZHAO Tiecheng,SHEN Yuefen,ZHU Guozhen. Temperature field calculation of boiler drum—finite element analysis of 3-D unsteady variable character uneven material heat conduction problem[J]. Proceedings of the CSEE,1997,17(4):217-220.
[9]CEBULA A,TALER J. Determination of transient temperature and heat flux on the surface of a reactor control rod based on temperature measurements at the interior points[J]. Applied Thermal Engineering,2014,63(1):158-169.
[10]JAREMKIEWICZ M,TALER D,SOBOTA T. Measuring transient temperature of the medium in power engineering machines and installations[J]. Applied Thermal Engineering,2009,29(16):3374-3379.
[11]MURIO D A.The mollification method and the numerical solution of the inverse heat conduction problem by finite differences[J].Computers & Mathematics with Applications,1989,17(10):1385-1396.
[12]TALER J,ZIMA W.Solution of inverse heat conduction problems using control volume approach[J].International Journal of Heat and Mass Transfer,1999,42(6):1123-1140.
[14]田廣,趙海飛.增壓鍋爐鍋筒內壓應力三維有限元分析[J].鍋爐制造,2013 (1):61-64.
TIAN Guang,ZHAO Haifei.Mechanical stress analysis of supercharged boiler drum by 3-D finite element method[J].Boiler Manufacturing,2013 (1):61-64.
A Coupling Method of Direct and Inverse Heat Conduction Problems for Transient Temperature Calculation of a Boiler Drum
LIBin,CHENFeng,SHILiangxiao
(School of Energy, Power and Mechanical Engineering, North China Electric Power University,Baoding 071003, Hebei Province, China)
A coupling method of direct and inverse heat conduction problems was proposed for transient temperature calculation of a pressurized boiler drum, of which the locally heated outer surface is divided into two regions, namely heated and unheated area. For the unheated area, the inverse method is used to determine the temperature field by circumferentially arranging thermocouples to measure the outside surface temperature of the drum, and for the heated area, the direct method is adopted to solve the temperature field; whereas at the junction of above two areas, the temperatures determined by the inverse method are used as known boundary conditions to calculate the temperature field of heated area by the direct method. By this way the coupling solution of direct and inverse heat conduction problems is realized, and the whole temperature distribution of the drum is obtained. A comparison was made to the drum temperature distribution during cold startup of the boiler respectively obtained by Ansys software and the coupling method. Results show that the coupling method of direct and inversion heat conduction problems has high precision and strong adaptability under complicated boundary conditions, which therefore may be used in actual engineering projects.
boiler drum; temperature field; inverse heat conduction problem; direct heat conduction problem; coupling method of direct and inverse heat conduction problem
A學科分類號:470.30
2014-02-17
2014-06-05
李斌(1969-),男,河北保定人,副教授,博士,主要從事傳熱過程數值模擬和熱力系統性能分析等方面的研究.
陳豐(通信作者),男,碩士研究生,電話(Tel.):15132461127;E-mail:chenfeng1513246@163.com.
1674-7607(2015)02-0096-07
TK223