鄭陸松趙潁杰孫寶芝齊洪亮
(1哈爾濱工程大學(xué)動力與能源工程學(xué)院,黑龍江 哈爾濱150001;2中國船舶重工集團公司第703研究所廣翰燃機事業(yè)部,黑龍江 哈爾濱150001)
在核電站實際運行過程中,伴隨著供電功率的變化,蒸汽發(fā)生器常需要降負(fù)荷運行[1]。運行經(jīng)驗表明,在降負(fù)荷過程中,一、二回路流體熱工水力特性呈現(xiàn)動態(tài)變化,并且易出現(xiàn)流動不穩(wěn)定性現(xiàn)象[2]。因此,針對蒸汽發(fā)生器降負(fù)荷運行過程中一、二次側(cè)流體流動特點,進行流體與傳熱管耦合傳熱的研究對蒸汽發(fā)生器安全有效運行有重要意義。
隨著計算機技術(shù)的發(fā)展,國內(nèi)外學(xué)者較多地采用數(shù)值模擬方法研究蒸汽發(fā)生器的熱工水力特性。早期的數(shù)值研究中,大多數(shù)學(xué)者采用多孔介質(zhì)模型模擬二次側(cè)汽、液兩相流沸騰傳熱過程[3-6]。多孔介質(zhì)模型沒有考慮一次側(cè)流體及傳熱管,不能準(zhǔn)確模擬一、二次側(cè)流體的耦合傳熱過程。針對多孔介質(zhì)模型的不足,相關(guān)學(xué)者基于已有模型進行改進[7-10]或者采用自主開發(fā)模型[11-12]計算蒸汽發(fā)生器熱工水力特性。然而,公開文獻中蒸汽發(fā)生器熱工水力特性的研究所采用的數(shù)值計算模型簡單,不能全面考慮一、二次側(cè)流動換熱過程,而且多為蒸汽發(fā)生器的穩(wěn)態(tài)數(shù)值研究,對變負(fù)荷過程中蒸汽發(fā)生器流動傳熱的動態(tài)特性研究較少。
本文以大亞灣核電站蒸汽發(fā)生器為原型,針對蒸汽發(fā)生器降負(fù)荷過程中流體流動傳熱特點,基于CFD技術(shù),采用CFX軟件對蒸汽發(fā)生器進行100%~50%工況的降負(fù)荷瞬態(tài)模擬。提出了一種能夠完整地計算蒸汽發(fā)生器內(nèi)一、二次側(cè)流體耦合流動傳熱特性的瞬態(tài)數(shù)值模擬方法。模擬結(jié)果揭示了蒸汽發(fā)生器流動傳熱特性在降負(fù)荷過程中的動態(tài)變化規(guī)律,所得結(jié)論可為蒸汽發(fā)生器安全運行提供參考。
基于實際蒸汽發(fā)生器,建立簡化的物理模型,示于圖1。模型包括一、二次側(cè)流體域,傳熱管及支撐板。U形傳熱管外徑為19.05mm,壁厚為1.09mm,管中心間距為27.43mm,直管高度為9m,彎管中心曲率半徑為0.82m,其幾何尺寸均與實際蒸汽發(fā)生器相同。一次側(cè)流體由U形傳熱管一端進入,沿管內(nèi)向上流動,經(jīng)彎管段后由傳熱管另一端向下流出。模型中一次側(cè)流體進入端為熱端,流出端為冷端。二次側(cè)流體由傳熱管兩端底部進入,沿管束間通道向上流動,由彎管段外側(cè)管束間流出。9塊四葉梅花形支撐板位于二次側(cè)流體域且沿傳熱管高度方向均勻分布以支撐傳熱管,二次側(cè)流體通過支撐板與傳熱管間的四葉流水孔流過支撐板,支撐板厚度為30mm。
根據(jù)所建立的物理模型,結(jié)合大亞灣核電站蒸汽發(fā)生器運行參數(shù)[13]設(shè)置100%工況的邊界條件。一回路進口質(zhì)量流量為4.16kg·s-1,進口溫度為327.6℃,出口平均靜壓15.5MPa。二回路工質(zhì)進口質(zhì)量流量為1.887kg·s-1,液相水溫度為270.3℃,出口平均靜壓為6.81MPa,并考慮重力影響。其中,二次側(cè)流體進口流量平均分配到兩個進口端,液相體積分?jǐn)?shù)為1,汽相體積分?jǐn)?shù)為0。U形傳熱管及支撐板壁面設(shè)置為流-固交界面,以模擬流體與壁面的耦合傳熱過程。采用變尺度壁面函數(shù)模型[14]模擬壁面區(qū)域湍流及水力結(jié)構(gòu),汽、水相分別為有、無滑移邊界條件。其余剖切面邊界設(shè)置為對稱面。
蒸汽發(fā)生器降負(fù)荷過程中,一、二次側(cè)流體邊界參數(shù)變化示于表1,其中的流體參數(shù)變化規(guī)律是根據(jù)實際測量參數(shù)[13]以及核電站運行標(biāo)準(zhǔn)規(guī)定的變化速率范圍[15]進行擬合?;趯嶋H運行過程,模擬中蒸汽發(fā)生器一次側(cè)進口流量及出口壓力在降負(fù)荷過程中保持不變。數(shù)值模擬中,首先根據(jù)流體邊界條件計算蒸汽發(fā)生器100%工況下的穩(wěn)態(tài)流動傳熱,然后再以100%工況的模擬結(jié)果為初始流場,并依據(jù)表1添加流體邊界參數(shù)的變化,進而模擬蒸汽發(fā)生器100%~50%工況的降負(fù)荷過程。CFX中總持續(xù)時間為4500s,其中,0s對應(yīng)100%工況,2700s對應(yīng)70%工況,4500s對應(yīng)50%工況。鑒于流體參數(shù)變化緩慢,一次側(cè)流體溫度變化速率僅為0.00392℃·s-1,所以采用時間步為0.1s,每時間步迭代50次,計算結(jié)果每100s保存一次。
圖1 蒸汽發(fā)生器簡化物理模型Fig.1 Simplified physical model of steam generator
表1 流體邊界參數(shù)變化擬合關(guān)系式Table 1 Fluid boundary parameters fitting equation
蒸汽發(fā)生器一次側(cè)為單相強制對流,由基本微分方程[16]描述,并采用標(biāo)準(zhǔn)k-ε模型描述湍流脈動。傳熱管的導(dǎo)熱過程由導(dǎo)熱微分方程控制。二次側(cè)流體為汽、液兩相流動且涉及沸騰相變傳熱,采用兩流體模型描述汽、液兩相流動與傳熱。兩流體模型中考慮了汽泡表面張力,采用Particle模型[17]計算兩相界面?zhèn)鬟f;采用熱相變模型[18]計算沸騰相變傳熱;液相湍流脈動采用標(biāo)準(zhǔn)k-ε模型,汽相采用離散相零方程模型,并采用STAO增強渦流黏度模型[19]描述兩相間湍流傳遞。由于一、二次側(cè)耦合流動的復(fù)雜性,相關(guān)理論模型可參考相應(yīng)文獻,而論文中不再展示。
圖2 簡化物理模型網(wǎng)格Fig.2 Grids of simplified physical model
蒸汽發(fā)生器簡化模型網(wǎng)格示于圖2。采用ICEM劃分流體計算網(wǎng)格,網(wǎng)格類型均為六面體結(jié)構(gòu)化網(wǎng)格。其中,一、二次側(cè)流體網(wǎng)格均在近壁處加密,以捕捉流動邊界。傳熱管沿壁厚方向分布3層網(wǎng)格,以計算一次側(cè)流體向二次側(cè)流體的耦合熱傳導(dǎo)。二次側(cè)流體網(wǎng)格在支撐板位置沿流動方向加密,以捕捉由于流動截面變化引起的湍流渦流流動。二次側(cè)流體在支撐板與傳熱管的縫隙位置采用3層網(wǎng)格加密,以捕捉縫隙區(qū)二次側(cè)流體的流動結(jié)構(gòu)。利用二次側(cè)流體出口平均質(zhì)量含汽率驗證網(wǎng)格無關(guān)解,示于圖3。當(dāng)流體網(wǎng)格數(shù)為4220000時達到無關(guān)解要求,相應(yīng)一次側(cè)流體域網(wǎng)格數(shù)為417792,二次側(cè)流體域網(wǎng)格數(shù)為3217980,傳熱管網(wǎng)格數(shù)為486048,支撐板網(wǎng)格數(shù)為98180。
圖3 不同網(wǎng)格數(shù)下二次側(cè)出口質(zhì)量含汽率Fig.3 Average steam quality of secondary outlet at different grids
根據(jù)蒸汽發(fā)生器100%~50%工況的降負(fù)荷瞬態(tài)模擬結(jié)果,對流體溫度、流速等熱工參數(shù)進行分析,初步探索了蒸汽發(fā)生器內(nèi)流動傳熱的動態(tài)特性。
一、二次側(cè)流體截面平均溫度沿管長方向的分布示于圖4。一次側(cè)流體溫度沿流動方向逐漸降低,且隨著工作負(fù)荷的減小,流體溫度降低的速率減小。對應(yīng)于100%、70%及50%工況,一次 側(cè) 流 體 出 口 溫 度 分 別 為 297.2、295.6 和294.8℃。計算結(jié)果與大亞灣核電站實際測量值基本相符[13],證明了數(shù)值模擬的準(zhǔn)確性。二次側(cè)流體溫度沿流動方向,在冷、熱端的預(yù)熱段均逐漸升高,在沸騰段流體溫度均保持在不同工況的飽和溫度附近不變。隨著工作負(fù)荷的降低,二次側(cè)流體飽和溫度升高。在預(yù)熱段的支撐板處,由于二次側(cè)流體流速加快,流體內(nèi)能向動能轉(zhuǎn)化,流體溫度有所降低。
如圖5所示,一次側(cè)流體截面溫度分布均勻,近壁處流體溫度由于傳熱而降低且存在溫度邊界層。二次側(cè)近壁流體溫度較高,遠(yuǎn)離傳熱管壁面的主流區(qū)溫度降低。在支撐板上側(cè)由于二次側(cè)流體湍流度增強,流體充分混合,強化了換熱,所以二次側(cè)流體截面溫度分布較支撐板下側(cè)均勻。支撐板受二次側(cè)流體傳熱影響,近壁區(qū)其表面溫度較高,遠(yuǎn)離傳熱管的區(qū)域其溫度與二次側(cè)流體溫度接近。
圖4 一、二次側(cè)流體溫度變化曲線Fig.4 Variation curves of fluid temperature in primary and secondary sides
圖5 一、二次側(cè)流體在支撐板處局部溫度分布云圖Fig.5 Local temperature distributions of primary and secondary side fluids at supporting plates
降負(fù)荷過程中,一次側(cè)流速隨空間及時間的變化示于圖6。沿管長方向,流體從進口位置到穩(wěn)定流動過程中,由于存在進口效應(yīng),流速有所增加,隨后由于換熱的進行流體溫度減小,密度增加,所以流速又逐漸減小。彎管段,一次側(cè)流體由于流動方向發(fā)生變化,并且在高流速下存在二次環(huán)流,所以流體流速在進入彎管段時快速減小,并且在流出彎管段時流速又快速回升。沿管長方向,流體流速變化速率隨著負(fù)荷的下降而減小,而且不同空間位置處流速隨時間的變化也逐漸變緩。
降負(fù)荷過程中,對應(yīng)于100%、70%及50%工況的二次側(cè)液相流速示于圖7。隨著負(fù)荷降低,二次側(cè)液相流速空間分布整體下降。二次側(cè)單相流區(qū),流體流速變化取決于進口流量,變化緩慢。在兩相流區(qū),隨著含汽率的不斷增加,冷、熱端流體流速沿高度方向均快速上升。支撐板位置處,由于流通截面減小,流體流速快速增加,并且隨著負(fù)荷降低,同一位置的流速峰值逐漸下降。單相流區(qū),支撐板對流體的加速作用隨負(fù)荷降低基本不變,而兩相流區(qū),隨著負(fù)荷降低,支撐板對流體的加速作用逐漸減弱。
圖6 降負(fù)荷運行過程一次側(cè)流體速度場分布Fig.6 Primary side fluid velocity field distribution in load drop running process
圖7 不同工況二次側(cè)液相流速沿管長方向變化曲線Fig.7 Secondary side fluid velocity curves along tube length direction under different conditions
支撐板處,二次側(cè)液相流速云圖隨時間的變化示于圖8。支撐板下側(cè)流體流速分布均勻,隨著負(fù)荷降低,其流速逐漸減小。流體流入支撐板流水孔后,隨著流通截面減小,流體相互擠壓,流速快速升高。t=0s,支撐板流水孔主流速度達到10m·s-1,且在支撐板上側(cè)形成射流。由于射流速度的差異,支撐板上側(cè)形成強烈的湍流渦流區(qū)域,而且隨著負(fù)荷降低,流水孔主流射流速度快速減小。如圖所示,t=1500s時,主流速度約為8m·s-1;t=3000s時,主流速度約為6.5m·s-1;而當(dāng)負(fù)荷降低到50%工況,即t=4500s時,主流速度減小到約4m·s-1。由圖中還可以看到,隨著負(fù)荷的降低,支撐板上側(cè)射流長度及渦流區(qū)范圍也在逐漸縮小。
圖8 支撐板處二次側(cè)流體流速瞬態(tài)云圖Fig.8 Secondary side fluid velocity counter at support plates
圖9 二次側(cè)質(zhì)量含汽率沿管長方向動態(tài)變化規(guī)律Fig.9 Dynamic changes of secondary side steam quality along tube length direction
二次側(cè)流體質(zhì)量含汽率隨時間及空間變化示于圖9。對應(yīng)于100%、70%及50%三種工況,二次側(cè)出口平均質(zhì)量含汽率分別為24.5%、16.7%和12.1%,相應(yīng)循環(huán)倍率為4~8,與大亞灣核電站實際運行過程的循環(huán)倍率相符[13]。支撐板位置處,由于流通截面減小,二次側(cè)流體流速升高,換熱加強,質(zhì)量含汽率在支撐板上側(cè)渦流區(qū)突然上升,并且隨著流動逐漸穩(wěn)定而又快速降低。隨著工作負(fù)荷降低,平均質(zhì)量含汽率整體值及沿高度方向的增長速率均減小,而且二次側(cè)冷、熱端的汽化位置均逐漸升高。t=0s,熱端汽化位置在第1塊支撐板處,冷端則在第3塊支撐板處。t=4500s,熱端汽化位置提高到了第2塊支撐板的位置,而冷端則提高到了接近第5塊支撐板的位置。二次側(cè)流體汽化后,含汽率沿流動方向均快速升高。冷、熱端流體進入彎管區(qū)域,含汽率均是先增加,而后伴隨著冷、熱端流體相遇并流出管束間換熱區(qū),含汽率突然降低。
支撐板處,二次側(cè)質(zhì)量含汽率分布云圖的瞬態(tài)變化示于圖10。隨著負(fù)荷降低,支撐板處局部截面的整體含汽率下降,而且渦流的影響范圍明顯縮小。在t=0s,即100%工況時,支撐板渦流區(qū)域局部含汽率高達0.67,而隨著負(fù)荷降低,在t=4500s時,相應(yīng)區(qū)域含汽率減小至不足0.3。由此可見,在蒸汽發(fā)生器降負(fù)荷過程中,支撐板處二次側(cè)流體相分布有較大改變。
圖10 支撐板處質(zhì)量含汽率瞬態(tài)變化云圖Fig.10 Transient variation counter of steam quality at support plates
不同工況下二次側(cè)流體表面換熱系數(shù)的變化曲線示于圖11。隨工作負(fù)荷降低,表面換熱系數(shù)的整體平均值減小,對應(yīng)于100%、70%及50%工 況, 其 平 均 值 分 別 為 45176.9、36319.3 和32458.9W·m-2·K-1,模擬結(jié)果符合Rohsenow經(jīng)驗公式[20]理論值,驗證了數(shù)值結(jié)果的準(zhǔn)確性。二次側(cè)冷、熱端表面換熱系數(shù)沿流動方向均呈遞增趨勢。單相對流區(qū),不同工況表面換熱系數(shù)基本相同,進入沸騰換熱區(qū),隨著含汽率的增加,表面換熱系數(shù)均快速增加。熱端傳熱強烈,沿流動方向表面換熱系數(shù)整體大于冷端。支撐板位置處,二次側(cè)流體表面換熱系數(shù)急劇增大又快速減小。冷、熱端流體在進入彎管區(qū)域時,流體表面換熱系數(shù)先增加,而在冷、熱端流體在彎管區(qū)相遇并流出管束區(qū)時,又快速減小。
圖12給出了支撐板處傳熱管外壁面換熱系數(shù)的瞬態(tài)變化。由于流水孔主流速度較高及射流渦流的影響,相應(yīng)區(qū)域表面換熱系數(shù)較大,且圍繞傳熱管呈周期分布。隨著負(fù)荷降低,相應(yīng)部位傳熱管的表面換熱系數(shù)逐漸減小。由圖8及圖12對比可知,表面換熱系數(shù)主要受二次側(cè)流體流速影響。
圖11 二次側(cè)表面換熱系數(shù)變化曲線Fig.11 Variation curves of surface heat transfer coefficient in secondary side
圖12 支撐板處傳熱管外壁表面換熱系數(shù)的瞬態(tài)變化Fig.12 Transient change counter of secondary side heat transfer coefficient at support plates
本文基于核電站實際運行過程,采用CFD技術(shù),對蒸汽發(fā)生器100%~50%工況降負(fù)荷過程進行瞬態(tài)模擬,揭示了蒸汽發(fā)生器流動傳熱特性的動態(tài)變化規(guī)律。所得主要結(jié)論如下:
(1)對應(yīng)于100%、70%及50%工況,一次側(cè)流體出口溫度分別為297.2、295.6和294.8℃,二次側(cè)出口平均質(zhì)量含汽率分別為24.5%、16.7%和12.1%,計算結(jié)果與大亞灣核電站實際值基本相符,而且二次側(cè)平均表面換熱系數(shù)分別為45176.9、36319.3和32458.9W·m-2·K-1,模擬結(jié)果符合Rohsenow經(jīng)驗公式的理論計算值。
(2)一次側(cè)流體溫度及流速沿流動方向均逐漸降低,而且隨著工作負(fù)荷的降低,相應(yīng)的變化速率均減小。沿流動方向,二次側(cè)流體溫度在預(yù)熱段逐漸升高,在沸騰段保持為飽和溫度,二次側(cè)流體流速在支撐板位置處突然增大又快速減小,且在支撐板上側(cè)出現(xiàn)渦流。隨著負(fù)荷的降低,二次側(cè)流體飽和溫度升高,流體主流速度減小,支撐板上側(cè)射流長度及渦流區(qū)范圍也在逐漸縮小。
(3)二次側(cè)流體質(zhì)量含汽率沿高度方向迅速升高,支撐板位置處,質(zhì)量含汽率突然上升又快速降低。降負(fù)荷過程中,二次側(cè)流體汽化位置延后,整體汽化強度降低,而且支撐板對二次側(cè)相分布的擾動作用也逐漸減弱。二次側(cè)流體表面換熱系數(shù)主要受二次側(cè)流體流速影響,沿流動方向呈遞增趨勢,且在支撐板處急劇增大又快速減小。隨著負(fù)荷降低,表面換熱系數(shù)減小。
[1] Li Haijun (李海軍).Operating characteristics analysis for steam generator of Ling Ao nuclear power station [D].Guangdong:South China University of Technology,2012:25-30.
[2] Hao Chengming (郝承明),Liu Xinkai(劉新凱),Peng Minjun(彭敏?。琗ia Genglei(夏庚磊),Liu Jiange (劉建閣).Research on influence to flow instability of OTSG under group operation mode [J].AtomicEnergyScienceand Technology(原子能科學(xué)技術(shù)),2011,45 (12):1466-1470.
[3] Jiang Xing (蔣興),Zhang Ming (張明),Xie Yongcheng(謝 永 誠 ),Yao Wenda (姚 偉 達 ).Three-dimensional numerical simulation of secondary side flow field in steam generator [J].AtomicEnergyScienceandTechnology(原子能科學(xué)技術(shù)),2008,42 (S1):438-443.
[4] Zhang Ronghua(張榮華),Xu Jijian(徐濟鑒),Chen Jun(陳軍).The thermal-h(huán)ydraulic simulation of compact reactor system [J].ChineseJournalofNuclearScienceand Engineering(核科學(xué)與工程),2004,24 (2):129-133.
[5] Ferng Y M,Chang H J.CFD investigating the impacts of changing operating conditions on the thermal-h(huán)ydraulic characteristics in a steam generator [J].AppliedThermal Engineering,2008,28:414-422.
[6] Ferng Y M,Ma Y P.Thermal-h(huán)ydraulic simulation of localized flow characteristics in a steam generator [J].Nuclear Technology,2000,136:186-194.
[7] Yang Zhilin (楊志林),Li Shixin(李師新),Xu Ming (徐明).The analysis of thermal hydraulies in heater of steam generator secondary side [J].JournalofShanghaiJiaoTongUniversity(上海交通大學(xué)學(xué)報),1996,30 (11):165-168.
[8] Hoeld A.A universally applicable thermal-h(huán)ydraulic driftflux based mixture-fluid 1Dmodel and code [J].Nuclear EngineeringandDesign,2007,237:1952-1967.
[9] Hoeld A.Application of the mixture-fluid channel element CCM within the U-tube steam generator code UTSG-3for transient thermal-h(huán)ydraulic calculations//15th International Conference on Nuclear Engineering [C].Nagoya,Japan,2007:22-26.
[10] Belliard M,Grandotto M.Computation of two-phase flow in steam generator using domain decomposition and local zoom methods[J].NuclearEngineeringandDesign,2002,213(2):223-239.
[11] Qiu Suizheng (秋 穗 正),Guo Yujun (郭 玉 君),Zhang Jinling(張金玲),Su Guanghui(蘇光輝),Yu Zhenwan(喻 真 烷 ),Jia Dounan (賈 斗 南 ).Dynamic simulating calculation for steam generator in nuclear power plant [J].JournalofNuclearScienceandEngineering(核科學(xué)與工程),1995,15 (4):289-295.
[12] Wang Chuan (王川),Yu Lei(于雷).Investigation on single phase water reverse flow in inverted U-tubes of steam generator under condition of natural circulation [J].Nuclear PowerEngineering(核動力工程),2011,32 (1):58-62.
[13] Ding Xunshen (丁訓(xùn)慎).Natural circulation steam generator operating principles and operating parameters [J].Nucl.PowerPlant(核電站),2003,(2):15-21.
[14] Yang Yuanlong (楊元龍).Numerical investigation on thermal-h(huán)ydraulic of steam generator based on two-fluid model[D].Harbin:Harbin Engineering University,2013:22-25.
[15] Xue Hanjun (薛漢?。?Nuclear Power Plant(核能動力裝置)[M].Beijing:Atomic Energy Press,2002:65-99.
[16] Zhang Shishuai (張 師 帥 ).Computational Fluid Dynamics and Application(計算流體動力學(xué)及其應(yīng)用)[M].Wuhan:Huazhong University of Science and Technology Press,2011:6-10,14-41.
[17] Wang Jiaqiong(汪家瓊),Jiang Wanming(蔣萬明),Kong Fanyu(孔繁余),Qu Xiaoyun(屈曉云),Su Xianghui(宿向 輝 ).Numerical simulation of solid-liquid two-phase turbulent flow in centrifugal pump based on particle model[J].JournalofDrainageandIrrigationMachinery Engineering(排 灌 機 械 工 程 學(xué) 報 ),2013,31 (10):846-850.
[18] Chen Erfeng,Li Yanzhong,Cheng Xianghua.CFD simulation of upward subcooled boiling flow of Refrigerant-113 using the two-fluid mode [J].AppliedThermalEngineering,2009,29:2508-2517.
[19] Sato Y,Sekoguchi K.Liquid velocity distribution in twophase bubble flow [J].Int.J.Multiphase.Flow,1975,2(1):79-95.
[20] Bjorg R W,Hall G R,Rohsenow W M.Correlation of forced convection boiling heat transfer data [J].InternationalJournalofHeatandMassTransfer,1982,25 (6):753-757.