張 淳，王富強(qiáng)，譚建宇，來慶志

(1．中國(guó)機(jī)械設(shè)備工程股份有限公司，北京 100055;

2．哈爾濱工業(yè)大學(xué)(威海)汽車工程學(xué)院，山東 威海 264209)

0 前言

國(guó)際上太陽能熱動(dòng)力發(fā)電的聚光系統(tǒng)主要有三種形式:槽式，塔式和碟式［1］。槽式聚光系統(tǒng)太陽能熱動(dòng)力發(fā)電技術(shù)的優(yōu)點(diǎn)是:容量可大可小;吸熱器和聚光系統(tǒng)都布置于地面上，安裝和維護(hù)方便;各聚光系統(tǒng)可同步跟蹤，控制系統(tǒng)成本低［2］。目前只有槽式聚光系統(tǒng)太陽能熱動(dòng)力發(fā)電技術(shù)已商業(yè)化運(yùn)行，而含有塔式聚光和碟式聚光太陽能熱動(dòng)力發(fā)電技術(shù)仍處于實(shí)驗(yàn)和示范階段［3］。

入射太陽光經(jīng)過槽式聚光系統(tǒng)反射后匯聚在管式吸熱器外表面的下半段。管式吸熱器外表面的下半段受到高匯聚、非均勻熱流密度的照射，而管式吸熱器外表面的上半段卻受到非匯聚的低太陽熱流密度的照射［4］。因此，管式吸熱器外表面的熱流密度場(chǎng)分布為高度非均勻的，管式吸熱器容易承受大溫度梯度并引起高熱應(yīng)力而導(dǎo)致吸熱器的失效［5］。

為了降低吸熱器的熱應(yīng)力，Verlotski 等通過控制流體速率達(dá)到降低吸熱器熱應(yīng)力的目的［6］;Flore等提出了一種銅－不銹鋼雙層管式吸熱器［7］;Lata等采用高鎳合金替代奧氏體不銹鋼作為管式吸熱器材料來實(shí)現(xiàn)熱應(yīng)力降低的目的［8］。采用實(shí)驗(yàn)方法研究太陽能吸熱器的熱應(yīng)力及熱變形是多數(shù)學(xué)者采用的方法;即便采用數(shù)值模擬的方法研究吸熱器的溫度場(chǎng)和熱應(yīng)力場(chǎng)的過程中，多數(shù)學(xué)者也是采用了簡(jiǎn)化的熱流密度場(chǎng)或溫度場(chǎng)作為熱分析的邊界條件［9］。為此，本文作者曾采用蒙特卡洛與有限元相耦合法分析了管式吸熱器的溫度場(chǎng)和熱應(yīng)力場(chǎng)［10］。

本文采用光熱力順序解耦計(jì)算法研究流體流速對(duì)管式吸熱器的溫度場(chǎng)、熱應(yīng)變以及熱應(yīng)力場(chǎng)的影響。相關(guān)計(jì)算結(jié)果可以為實(shí)際運(yùn)行過程中防止運(yùn)行溫度過高以及對(duì)熱應(yīng)變和熱應(yīng)力的抑制提供參考。

1 光熱力順序解耦計(jì)算法

在光熱力順序解耦計(jì)算過程中，管式吸熱器外表面熱流密度場(chǎng)的計(jì)算采用本文作者等根據(jù)蒙特卡洛法編寫的程序代碼進(jìn)行計(jì)算［11］。蒙特卡洛法作為一種概率統(tǒng)計(jì)的數(shù)值計(jì)算方法，具有物理概念清晰、對(duì)復(fù)雜問題適應(yīng)性強(qiáng)等優(yōu)點(diǎn)而被計(jì)算物理學(xué)界廣泛應(yīng)用。采用蒙特卡洛法求解太陽能聚集傳輸問題的基本思想是:將太陽能看作是大量的、獨(dú)立的太陽光線組成;為了保證太陽光線分布的均勻性，假設(shè)每根光線攜帶相同的能量;將每根光線的傳輸過程分解為發(fā)射、反射、透射及吸收等一系列相互獨(dú)立的子過程，每個(gè)子過程都遵循特定的分布函數(shù)概率模型［12－13］。

在整個(gè)熱流密度場(chǎng)計(jì)算過程中，采用自編網(wǎng)格節(jié)點(diǎn)劃分程序?qū)⒐苁轿鼰崞鞯耐獗砻骐x散成軸向?yàn)?00 個(gè)節(jié)點(diǎn)、徑向?yàn)?00 個(gè)節(jié)點(diǎn)的網(wǎng)格。在光熱力順序解耦計(jì)算過程中熱流密度場(chǎng)計(jì)算的邊界條件為:太陽光不平行度為16';太陽輻照度為1 000 W/m2;環(huán)日比為0．05;忽略跟蹤誤差、指向誤差及位置誤差等對(duì)熱流密度場(chǎng)分布的影響。在光熱力順序解耦計(jì)算過程中采用的槽式聚光系統(tǒng)、管式吸熱器及玻璃罩的幾何尺寸和物性參數(shù)如表1 所示。在實(shí)際應(yīng)用中，多根管式吸熱器串聯(lián)到一起可達(dá)上百米。為了計(jì)算簡(jiǎn)化，本文僅取流體入口處2 m 長(zhǎng)管式吸熱器作為研究對(duì)象。

表1 槽式聚光系統(tǒng)和管式吸熱器的幾何尺寸及物性參數(shù)

為了增加管式吸熱器的玻璃罩對(duì)入射太陽能光線的透射率以及降低吸熱器熱損失，管式吸熱器的玻璃罩會(huì)敷有選擇性透過涂層。此外，管式吸熱器的金屬表面還會(huì)涂以選擇性吸收涂層。管式吸熱器金屬表面的選擇性吸收涂層在太陽光譜下的高吸收率和工作溫度下的低發(fā)射率不但可以提高吸熱器對(duì)入射太陽光的吸收，還可以有效的降低吸熱器的輻射換熱損失。由表1 可知，管式吸熱器玻璃罩的透過率高達(dá)0．965，非常接近于1．0，而且玻璃罩的厚度非常薄，僅為1．9 ×10－3m［14］，玻璃罩對(duì)管式吸熱器金屬表面的熱流密度場(chǎng)的數(shù)值大小及分布影響非常小，因此本文在光熱力順序解耦計(jì)算過程中將不考慮玻璃罩對(duì)管式吸熱器表面熱流密度分布的影響。

本文采用蒙特卡洛法計(jì)算得到的管式吸熱器外表面下半段的熱流密度分布如圖1 所示。為了將計(jì)算得到的熱流密度導(dǎo)入到溫度場(chǎng)計(jì)算中，采用擬合函數(shù)法將計(jì)算得到的管式吸熱器外表面下半段的熱流密度分布進(jìn)行多段函數(shù)擬合，并作為管式吸熱器溫度場(chǎng)分析過程中的第二類邊界條件。擬合函數(shù)得到的熱流密度場(chǎng)與計(jì)算得到的熱流密度場(chǎng)吻合良好，擬合誤差小于0．01%［10］。

將溫度場(chǎng)分析得到的網(wǎng)格節(jié)點(diǎn)溫度場(chǎng)通過插值的形式映射到有限元熱應(yīng)力分析網(wǎng)格節(jié)點(diǎn)上，并作為體力載荷施加在后序的熱應(yīng)力場(chǎng)分析中，并作為熱應(yīng)力場(chǎng)分析的輸入載荷。本文對(duì)管式吸熱器的熱應(yīng)力場(chǎng)分析過程中無外部邊界約束，因此管式吸熱器熱應(yīng)力的產(chǎn)生是由于吸熱器壁面上不均勻溫度分布而引起的。關(guān)于光熱力順序解耦計(jì)算法的更詳細(xì)介紹及相關(guān)的模型驗(yàn)證，請(qǐng)參見文獻(xiàn)［5，10－11，13］。

圖1 管式吸熱器外表面下半段的熱流密度場(chǎng)分布

2 計(jì)算網(wǎng)格劃分

在采用有限元進(jìn)行熱應(yīng)力場(chǎng)分析的過程中，模型的網(wǎng)格劃分質(zhì)量對(duì)求解結(jié)果的精度有著重要的影響;提高網(wǎng)格劃分質(zhì)量，減少離散化隨意性帶來的誤差是有限元分析過程中的重要工作。原則上，Ansys有限元分析只要求網(wǎng)格劃分質(zhì)量在0．2 以上即可計(jì)算求解;但是有限元分析過程中的網(wǎng)格質(zhì)量不高會(huì)帶來計(jì)算時(shí)間長(zhǎng)、收斂穩(wěn)定性差、甚至計(jì)算結(jié)果不準(zhǔn)確等缺點(diǎn)［10］。

圖2 溫度場(chǎng)分析和熱應(yīng)力場(chǎng)分析的網(wǎng)格劃分

為了更加精確的捕捉熱應(yīng)力沿吸熱器管壁徑向和切向方向的分布，在熱應(yīng)力分析過程中的管壁網(wǎng)格劃分密度比溫度場(chǎng)分析過程中的管壁網(wǎng)格劃分密度高。本文采用Ansys Workbench 中提供的ICEM網(wǎng)格劃分模塊對(duì)管式吸熱器的流體及管壁幾何模型進(jìn)行網(wǎng)格劃分，網(wǎng)格劃分過程中采用映射網(wǎng)格技術(shù)產(chǎn)生O 型結(jié)構(gòu)化四面體網(wǎng)格。網(wǎng)格劃分結(jié)果如圖2所示，其中用于溫度場(chǎng)計(jì)算分析的網(wǎng)格劃分單元數(shù)為95 888，用于熱應(yīng)力場(chǎng)網(wǎng)格劃分單元數(shù)為175 168。在進(jìn)行溫度場(chǎng)和應(yīng)力場(chǎng)分析前的網(wǎng)格無關(guān)性檢查中發(fā)現(xiàn)，本文劃分的網(wǎng)格單元數(shù)足以保證計(jì)算精度。網(wǎng)格質(zhì)量檢查結(jié)果表明用于溫度場(chǎng)分析的網(wǎng)格劃分質(zhì)量均在0．65 以上;用于熱應(yīng)力場(chǎng)分析的網(wǎng)格劃分質(zhì)量更高，均在0． 95 以上，遠(yuǎn)遠(yuǎn)高于Ansys 有限元分析的網(wǎng)格劃分質(zhì)量要求。

3 流速對(duì)溫度場(chǎng)分布的影響

圖3 所示為不同流速下管式吸熱器出口端面處溫度場(chǎng)分布。由圖3 可知，隨著流體流速的增加，管式吸熱器出口端面處的溫度迅速的降低;流體未被加熱區(qū)域也不斷的增大。當(dāng)流體流速為0． 1 m/s時(shí)，管式吸熱器出口端面處的最高溫度為430 K;而當(dāng)流體流速增加到0．5 m/s 時(shí)，管式吸熱器出口端面處的最高溫度僅為350 K。

圖4 所示為不同流速下管式吸熱器內(nèi)壁面溫度場(chǎng)沿長(zhǎng)度方向分布。由圖可知，對(duì)于流體流速為0.3 m/s 和0．5 m/s 兩種計(jì)算工況，入口段的溫度分布沿吸熱器長(zhǎng)度方向迅速地增加(z≤0．1 m)，而后溫度沿長(zhǎng)度方向的增加速度變緩，溫度分布趨于線性增加;對(duì)于流體流速為0．1 m/s 的計(jì)算工況，在整個(gè)管式吸熱器長(zhǎng)度下溫度場(chǎng)沿長(zhǎng)度方向都是迅速增加。

圖3 管式吸熱器出口端面處溫度場(chǎng)分布

4 流速對(duì)總應(yīng)變分布的影響

圖5 所示為不同流速下管式吸熱器出口端面處的總應(yīng)變沿圓周角分布。由圖可知，管式吸熱器的總應(yīng)變隨著流體流速的增加而迅速的降低，當(dāng)流體流速由0．1 m/s 增加到0．5 m/s 時(shí)，總應(yīng)變由4．14×10－4降低到0．62 ×10－4。在流體流速為0．1 m/s和0．3 m/s 兩種計(jì)算工況時(shí)，總應(yīng)變沿圓周角分布趨勢(shì)是一致的:總應(yīng)變的最大值都出現(xiàn)在圓周角為270°的位置處，兩條曲線都在圓周角為38°、140°、234°及305°四個(gè)位置處出現(xiàn)了曲線的拐點(diǎn)。流體流速為0．1 m/s 管式吸熱器的總應(yīng)變與流體流速為0．3 m/s 管式吸熱器的總應(yīng)變相差不大，流體流速為0．1 m/s 時(shí)的管式吸熱器總應(yīng)變比流體流速為0．3 m/s 的管式吸熱器總應(yīng)變稍高一點(diǎn)，總應(yīng)變最高差值為22%，出現(xiàn)在圓周角為78°位置處。

圖4 管式吸熱器內(nèi)壁面溫度沿長(zhǎng)度分布

圖5 管式吸熱器出口端面處總應(yīng)變沿圓周角分布

5 流速對(duì)熱應(yīng)力場(chǎng)的影響

圖6 所示為不同流體流速下管式吸熱器內(nèi)壁面軸向熱應(yīng)力沿長(zhǎng)度方向分布。由圖可知，軸向熱應(yīng)力隨著流體流速的增加而降低。在三種不同流體流速計(jì)算情況下，管式吸熱器入口端面處的軸向熱應(yīng)力均呈現(xiàn)為拉應(yīng)力且數(shù)值非常小，接近于零;隨著長(zhǎng)度的不斷增加，軸向熱應(yīng)力迅速的由拉應(yīng)力轉(zhuǎn)變?yōu)閴簯?yīng)力，并且達(dá)到一個(gè)很高的、接近各自曲線最大值的壓應(yīng)力;在管式吸熱器出口端面處，軸向熱應(yīng)力迅速的降低到接近于零。引起這種現(xiàn)象的原因如下:管式吸熱器外表面圓周角為270°的位置處受到的熱流密度最大，溫度也最高;管式吸熱器外表面圓周角為270°的位置受到高溫后的膨脹速度比周邊位置膨脹速度快，因此圓周角為270°的位置會(huì)受到來自周邊溫度相對(duì)低、膨脹速度相對(duì)慢的位置的約束，造成管式吸熱器中間段區(qū)間內(nèi)的軸向熱應(yīng)力為壓應(yīng)力而非拉應(yīng)力;而在管式吸熱器兩端面處，兩個(gè)端面為自由膨脹端，沒有外界約束，可以自由的沿軸向方向的兩端面處膨脹，因此兩端面處軸向熱應(yīng)力比較小，而且接近于零。

圖6 管式吸熱器內(nèi)壁面軸向熱應(yīng)力沿長(zhǎng)度分布

圖7 管式吸熱器內(nèi)壁面徑向熱應(yīng)力沿長(zhǎng)度分布

圖7 所示為不同流體流速下管式吸熱器徑向熱應(yīng)力沿長(zhǎng)度分布。由圖可知，徑向熱應(yīng)力隨著流體流速的增加而降低。三種不同流體流速計(jì)算情況時(shí)，管式吸熱器的徑向熱應(yīng)力隨著長(zhǎng)度的增加，迅速的由入口端面處的拉應(yīng)力變?yōu)閴簯?yīng)力;在管式吸熱器長(zhǎng)度為0．1 m 到1．9 m 區(qū)間內(nèi)，三種不同流體流速情況時(shí)管式吸熱器的徑向熱應(yīng)力均變化不大，維持在接近于零的壓應(yīng)力。在管式吸熱器出口端面處，三種流體流速情況時(shí)管式吸熱器的徑向熱應(yīng)力又迅速的由壓應(yīng)力變?yōu)槔瓚?yīng)力。三種不同流體流速情況時(shí)的管式吸熱器最大徑向熱應(yīng)力都出現(xiàn)在吸熱器入口處。

圖8 所示為不同流體流速下管式吸熱器切向熱應(yīng)力沿長(zhǎng)度分布。在管式吸熱器兩端面處的切向熱應(yīng)力的數(shù)值遠(yuǎn)高于沿吸熱器長(zhǎng)度方向其他位置處的切向熱應(yīng)力的數(shù)值。在管式吸熱器入口端面處，切向熱應(yīng)力表現(xiàn)為拉應(yīng)力，沿著吸熱器長(zhǎng)度的增加，切向熱應(yīng)力迅速的由拉應(yīng)力轉(zhuǎn)變?yōu)閴簯?yīng)力。在管式吸熱器長(zhǎng)度0．1 m 到1．9 m 的范圍內(nèi)，三種不同流體流速情況時(shí)管式吸熱器的切向熱應(yīng)力均表現(xiàn)為壓應(yīng)力，且三種不同流速情況時(shí)管式吸熱器的切向熱應(yīng)力都各自維持在一個(gè)穩(wěn)定值附近變化。在相同位置處，切向熱應(yīng)力隨著流體流速的增加而降低。與軸向熱應(yīng)力和徑向熱應(yīng)力相比，切向熱應(yīng)力的最大值大于軸向熱應(yīng)力和徑向熱應(yīng)力。

圖8 管式吸熱器內(nèi)壁面切向熱應(yīng)力沿長(zhǎng)度分布

圖9 管式吸熱器出口端面處等效熱應(yīng)力沿圓周角分布

由形變應(yīng)變能理論可知，等效熱應(yīng)力與徑向熱應(yīng)力、軸向熱應(yīng)力和切向熱應(yīng)力的關(guān)系為［15］

式中 σr、σz、σφ和σeff——徑向熱應(yīng)力、軸向熱應(yīng)力、切向熱應(yīng)力和等效熱應(yīng)力。

圖9 所示為不同流體流速下管式吸熱器的出口端面處等效熱應(yīng)力沿圓周角分布。由圖可知，管式吸熱器的等效熱應(yīng)力隨流體流速的增加而迅速降低。當(dāng)流體流速由0．1 m/s 增加到0．5 m/s 時(shí)，管式吸熱器的等效熱應(yīng)力由72． 7 MPa 降低到10．9 MPa。與總應(yīng)變沿圓周角分布規(guī)律一樣，在流體流速為0．1 m/s 和0．3 m/s 兩種計(jì)算工況下，等效熱應(yīng)力沿圓周角分布趨勢(shì)是一致的，等效熱應(yīng)力的最大值都出現(xiàn)在圓周角為270°的位置處，兩條曲線都在圓周角為38°、140°、234°及305°四個(gè)位置處出現(xiàn)了曲線的拐點(diǎn)。根據(jù)公式(1)可知，與軸向熱應(yīng)力和徑向熱應(yīng)力相比，切向熱應(yīng)力對(duì)管式吸熱器的等效熱應(yīng)力貢獻(xiàn)比例更大。因此，均勻化溫度沿管式吸熱器圓周方向的分布是熱應(yīng)力抑制的有效手段之一。

綜合以上分析，本文數(shù)值計(jì)算結(jié)果與Verlotski［6］的實(shí)驗(yàn)結(jié)論相類似，通過控制流體的流速可以有效的降低管式吸熱器的總應(yīng)變及等效熱應(yīng)力。

6 結(jié)論

本文采用光熱力順序解耦計(jì)算法分析了流體流速對(duì)管式吸熱器溫度場(chǎng)、總應(yīng)變及熱應(yīng)力場(chǎng)的影響。計(jì)算結(jié)果表明管式吸熱器的溫度場(chǎng)、總應(yīng)變及等效熱應(yīng)力均隨著流體流速的增加而降低，控制流體速率可以有效的降低管式吸熱器的總應(yīng)變和熱應(yīng)力。與軸向熱應(yīng)力和徑向熱應(yīng)力相比，切向熱應(yīng)力對(duì)管式吸熱器的等效熱應(yīng)力的貢獻(xiàn)比例更大。因此，均勻化溫度沿管式吸熱器的圓周方向分布是熱應(yīng)力抑制的有效手段之一。

［1］屠楠，魏進(jìn)家，方嘉賓．太陽能腔式吸熱器表面反射率的選擇［J］．工程熱物理學(xué)報(bào)，2014，35(4):700 －704．

［2］何梓年．太陽能熱利用［M］．合肥:中國(guó)科學(xué)技術(shù)大學(xué)出版社，2009．

［3］Mao Q J，Xie M，Tan H P． Effects of material selection on the radiation flux of tube receiver in a dish solar system［J］．Heat Transfer Research，2014，45(4):339 －347．

［4］Cheng Z D，He Y L，Xiao J，Cui F Q，Du B C，Zheng Z J，Xu Y． Comparative and sensitive analysis for parabolic trough solar collectors with a detailed Monte Carlo ray － tracing optical model［J］．Appl． Energy，2014(115):559－572．

［5］Wang F Q，Shuai Y，Yuan Y，Yang G，Tan H P．Thermal stress analysis of eccentric tube receiver using concentrated solar radiation［J］．Solar Energy，2010，85(10):1809－1815．

［6］Verlotski V，Schaus M，Pohl M． A solar thermal mgopowder receiver with working temperatures of more than 1600℃:model investigation by using laser as an irradiation source［J］．Solar Energy Materials and Solar Cells，1997，45(3):227－239．

［7］Flores V C，Almanza R F． Behavior of compound wall copper－steel receiver with stratified two－phase flow regimen in transient states when solar irradiance is arriving on one side of receiver［J］．Solar Energy，2004(76):195 －198．

［8］Lata J M，Rodriguez M A，Lara M A． High flux central receivers of molten salts for the new generation of commercial stand－alone solar power plants［J］． Journal of Solar Energy Engineering，2008(130):0211002．

［9］Kumar N S，Reddy K S． Thermal analysis of solar parabolic collector with porous disc receiver［J］． Applied Energy，2009(86):1804 －1812．

［10］Wang F Q，Shuai Y，Yuan Y，Liu B． Effects of material selection on the thermal stresses of tube receiver under concentrated solar irradiation［J］．Materials and Design，2012(33):284 －291．

［11］Wang F Q，Shuai Y，Tan H P，Lin R Y，Cheng P． Researches on a new type of solar surface cladding reactor with concentration quartz window［J］．Solar Energy，2013(94):177－181．

［12］Dai G L，Xia X L，Hou G F． Transmission performances of solar windows exposed to concentrated sunlight［J］．Solar Energy，2014(103):125 －133．

［13］Wang F Q，Tan J Y，Ma L X，Wang C C． Effects of glass cover on heat flux distribution for tube receiver with parabolic trough collector system［J］． Energy Conversion and Management，2015(90):47 －52．

［14］Almanza R F，Lentz A，Jimenez G． Receiver behavior in direct steam generation with parabolic toughs［J］． Solar Energy，1997(61):275 －278．