壓水堆蒸汽發(fā)生器三維兩相熱工水力分析程序開發(fā)與驗證

(上海核工程研究設計院有限公司，上海 200233)

0 引言

蒸汽發(fā)生器作為核電廠一、二回路的樞紐，是核島最重要的設備之一。其是否安全運行直接關系到整個核電廠的安全。蒸汽發(fā)生器二次側(cè)流場數(shù)據(jù)是蒸汽發(fā)生器設計和傳熱管流致振動和磨損分析的重要輸入。由于蒸汽發(fā)生器結(jié)構(gòu)復雜、體積龐大，采用試驗手段獲得二次側(cè)流場數(shù)據(jù)代價太大，難以實現(xiàn)；同樣由于結(jié)構(gòu)復雜，直接針對管束建模分析很難實現(xiàn)。因此本文將二次側(cè)管束區(qū)作為多孔介質(zhì)，引入分布阻力概念，對蒸汽發(fā)生器二次側(cè)三維兩相熱工水力分析程序進行開發(fā)，并開展小規(guī)模管束的流動換熱試驗，以完成對程序的試驗驗證。

最早將多孔介質(zhì)模型應用于換熱器數(shù)值模擬分析的是Patankar等[1]，他們于1974年將管殼式換熱器的殼側(cè)作為多孔介質(zhì)，完成了對管殼式換熱器殼側(cè)流場的數(shù)值模擬。在此基礎上，Sha 等[2]模擬了SG和反應堆堆芯中冷卻劑的流動；Prithiviraj等[3-4]模擬了三維換熱器中的流動。Ferng等[5-7]采用CFX 軟件多孔介質(zhì)模型對SG管束區(qū)的二次側(cè)流場進行了分析。1977年美國電力研究院(EPRI)開始開發(fā)蒸汽發(fā)生器熱工水力計算程序，后來通過試驗驗證和實堆反饋數(shù)據(jù)的不斷修正，于1984年正式推出了適用于U形管蒸汽發(fā)生器(UTSG)和直流蒸汽發(fā)生器(OTSG)二次側(cè)三維流場分析的CFD程序ATHOS[8]。Chan等[9]于1986年開發(fā)了UTSG穩(wěn)態(tài)、瞬態(tài)三維流場分析程序PORTHOS。法國電力公司(EDF)為研究PWR和LMFBR堆芯和組件內(nèi)的單相、兩相流場而開發(fā)了THYC程序[10]，之后推出適用于傳熱器流場分析的版本THYC-EXCHANGER[11]來分析蒸汽發(fā)生器二次側(cè)管束區(qū)兩相流場。法國原子能委員會(CEA)開發(fā)并驗證了UTSG三維熱工水力分析程序GENEPI[12]?？梢钥闯觯瑖H上已有多個蒸汽發(fā)生器三維熱工水力分析程序，這些軟件均基于多孔介質(zhì)模型，只是其兩相流控制方程、管內(nèi)外傳熱、相變及流動阻力計算方法有所不同。國內(nèi)也有學者[13-16]采用CFD軟件對蒸汽發(fā)生器二次側(cè)進行三維熱工水力分析，大多屬于方法性的探索研究，未形成專用分析程序。因此，本文基于已有的理論研究基礎，彌補現(xiàn)有分析程序的缺陷，開發(fā)可用于壓水堆蒸汽發(fā)生器三維兩相熱工水力分析程序。

1 程序開發(fā)

與常規(guī)流體的兩相流動類似，宏觀尺度的多孔介質(zhì)內(nèi)的兩相流動研究方法也分為三種：均相流模型、漂移流模型和兩流體模型。在早期的研究工作中，為簡化計算，多采用均相流模型，如CALIPSOS[17-18]、THEDA-1[19]、THEDA-2[20-21]等程序以及Ferng等[22-23]的一系列工作。均相流模型將兩相流體等效為單相流體，不能描述兩相間的相對速度，因此在計算豎直通道內(nèi)兩相流動的時候，會出現(xiàn)空泡份額偏大的問題。從理論上來說，兩流體模型可以準確地描述兩相間的熱力學和水力學不平衡特性，但由于管束外流動的相間作用力及相間傳熱等模型目前尚不完備，采用該模型需引入較多假設，導致兩流體模型的計算精度不能保證。漂移流模型同樣將兩相流體視為混合相，但是為描述相間的速度差，在均相方程的基礎上引入相間滑移速度。雖然漂移流模型從理論上相對于兩流體模型精度較差，但由于該模型需引入的假設條件較少，計算精度可滿足SG內(nèi)模擬的要求[8]。因此本程序開發(fā)中采用漂移流模型。多孔介質(zhì)內(nèi)的漂移流模型控制方程如下。

質(zhì)量守恒方程：

(1)

式中t——時間，s；

β——孔隙率；

ρm——混合物密度，kg/m3,ρm=αgρg+αlρl；

動量守恒方程：

(2)

式中p——壓力，Pa；

μm,eff——混合物有效黏度，Pa·s,

μm=αgμg+αlμl;

αg——汽相體積份額；

ρg,ρl——汽相、液相密度，kg/m3；

αl——液相體積份額；

能量守恒方程：

(3)

式中Hg——汽相焓，J/kg；

Hl——液相焓，J/kg；

km——混合物導熱系數(shù)，W/(m·K),km=αgkg+αlkl;

SE——能量源項,W/m3。

空泡份額輸運方程：

(4)

式中Sg——汽相質(zhì)量源項,kg/(m3·s)。

程序動量方程源項中考慮的阻力主要由以下部件引入，如傳熱管束、汽水分離器、下降通道、支承板等，其中管束、汽水分離器和下降通道的阻力以分布阻力形式添加到這些部件所在網(wǎng)格中，支承板阻力以集中阻力形式添加到支承板所在位置的網(wǎng)格界面上。根據(jù)流動方向不同，將流動方向沿管束軸向和橫向分解，在軸向和橫向分別引入動量源項。對于順流和橫掠阻力分別采用管束外順流和橫掠[24-25]阻力經(jīng)驗關系式計算。支承板、下降通道和汽水分離器的阻力系數(shù)，根據(jù)蒸汽發(fā)生器設計結(jié)構(gòu)直接給出。

能量方程的源項，即一次側(cè)向二次側(cè)釋熱，被簡化為一維分布計算。程序采取一維-三維耦合傳熱計算方式。該方法一、二次側(cè)網(wǎng)格示意如圖1所示，計算思路如下：將一次側(cè)沿流體流動方向劃分為J個控制體，局部網(wǎng)格編號用j表示；將二次側(cè)劃分為K個控制體，局部網(wǎng)格編號為k；通過DEFINE_INIT宏中的幾何處理模塊、根據(jù)網(wǎng)格位置將一、二次側(cè)網(wǎng)格進行匹配。

一次側(cè)為單相液體對流換熱，采用Dittus-Boelter 公式[26]計算；二次側(cè)換熱包括單相對流換熱、過冷沸騰換熱及飽和沸騰換熱，對于單相對流換熱，由于目前的經(jīng)驗關系式多為管外順流和橫掠關系式，無傾斜沖刷管束的換熱關系式，因此將流速分解為順流和橫掠流速，分別計算順流換熱系數(shù)[24-25]和橫掠換熱系數(shù)[27]，并將兩者的加權和[8]作為實際的斜掠換熱系數(shù)。在計算中，給定一次側(cè)入口的質(zhì)量流速和溫度，便可計算一、二次側(cè)間的流動換熱。

圖1 程序的一、二次側(cè)網(wǎng)格示意

考慮到多孔介質(zhì)模型對局部參數(shù)的平均效應，以及漂移流模型對相界面密度的依賴性較低，因此，程序中對沸騰模型進行簡化，僅根據(jù)能量平衡計算兩相流溫度及蒸發(fā)率，不計算相界面輸運過程。引入熱量密度二分方法，將總熱流密度分為蒸發(fā)熱流密度和對流熱流密度，分別用于蒸發(fā)液相和加熱混合相。

空泡份額方程中的源項表示氣相質(zhì)量源項，其值等于蒸發(fā)率與冷凝率之差，對于穩(wěn)態(tài)計算忽略SG內(nèi)的非穩(wěn)態(tài)冷凝現(xiàn)象，氣相質(zhì)量源項等于液相蒸發(fā)率。

對于湍流模型，已公開發(fā)表的文獻多是針對純流體區(qū)域的精細網(wǎng)格提出[28]，Teruel等[29-31]從宏觀尺度研究多孔介質(zhì)內(nèi)的湍流流動，他們在宏觀網(wǎng)格內(nèi)對k-e方程進行積分，得到宏觀k-e二方程湍流模型，并且分別采用微觀湍流模型和宏觀湍流模型計算幾何規(guī)則的多孔介質(zhì)內(nèi)的流動，從而驗證這些宏觀k-e湍流模型的正確性。但是這些模型均針對多孔介質(zhì)內(nèi)的單相流體進行開發(fā)，并不適用于多孔介質(zhì)內(nèi)的兩相流動。因此，本程序開發(fā)中采用了商用蒸汽發(fā)生器熱工水力分析程序ATHOS[32]和GENEPI[33]中采用的零方程湍流模型計算兩相混合物的有效黏度。

針對壓水堆蒸汽發(fā)生器結(jié)構(gòu)特點，以CFD軟件多孔介質(zhì)模型為基礎，采用UDF將以上傳熱、阻力模型添加到求解器中，并通過UDF求解一次側(cè)流場及一、二次側(cè)耦合傳熱特性，完成蒸汽發(fā)生器二次側(cè)三維兩相熱工水力分析程序的開發(fā)。

2 程序驗證

在本文開發(fā)的程序用于工程之前，需對其進行試驗驗證。由于蒸汽發(fā)生器結(jié)構(gòu)龐大，對原型蒸汽發(fā)生器進行試驗研究需要極大的蒸汽流量和加熱功率，成本很高。因此，基于?；治?，采用小規(guī)模管束進行流動傳熱試驗，模擬蒸汽發(fā)生器一、二次側(cè)流場，獲得試驗件一、二次側(cè)流場信息，同時采用該程序?qū)υ囼灱７治?，預測試驗工況下的試驗件的一、二次流場數(shù)據(jù)，通過與試驗測量結(jié)果的比較，達到對程序驗證的目的。

2.1 蒸汽發(fā)生器二次側(cè)流動傳熱試驗

試驗中采用40根傳熱管模擬二次側(cè)管束的流動換熱狀態(tài)，管束橫截面如圖2所示。為準確體現(xiàn)蒸汽發(fā)生器原型的熱工水力參數(shù)，試驗件管束高度上與原型保持相等，試驗模擬體中采用的傳熱管材料、管徑、管間距、布置方式、管束長度和支承板結(jié)構(gòu)與原型相同。試驗中一、二次側(cè)進口溫度、系統(tǒng)運行壓力和循環(huán)倍率等參數(shù)與原型相同，分別為：傳熱管規(guī)格為?17.48 mm×1.01 mm，直段高度9.2 m，相鄰傳熱管的中心距24.89 mm，呈三角形排列，一次側(cè)壓力15.5 MPa，進口溫度322 ℃；二次側(cè)壓力5.6 MPa，給水溫度226.7 ℃。試驗模擬體如圖3所示，包括上下筒體、上下封頭、管板、管束、下降管、管束支承板、汽水分離裝置等部件。為了避免測量探頭同時穿過筒體和套筒兩層壁面，采用兩根圓形下降管代替蒸汽發(fā)生器環(huán)形下降通道，同時在下降管上安裝了調(diào)節(jié)球閥，調(diào)節(jié)下降管阻力，以便更好地控制試驗模擬體運動的循環(huán)倍率，下降管與管束區(qū)入口的連接如圖4所示。試驗體運行及設計參數(shù)如表1,2所示。

圖2 管束橫截面示意

圖3 試驗模擬體示意

圖4 下降管與管束區(qū)入口連接示意

表1 試驗體額定運行參數(shù)

表2 試驗體設計參數(shù)

試驗中測量的參數(shù)包括流量、溫度、壓力、差壓、空泡份額、流體速度等。一次側(cè)進、出口溫度通過在一次側(cè)進、出口管道上安裝鎧裝熱電偶測量；一次側(cè)流量通過在一次側(cè)進口管道上安裝孔板流量計測量。二次側(cè)出口蒸汽壓力在蒸發(fā)器出口管上利用壓力傳感器測量，采用在蒸汽出口管道上安裝的孔板流量計測量蒸汽出口流量，給水流量測量與蒸汽流量測量方法一致，在給水管道上安裝鎧裝熱電偶測量給水溫度。采用光纖探針和伽馬射線法測量空泡份額。通過測量得到如下宏觀參數(shù)：熱功率、一回路壓力、一回路進、出口溫度、給水流量、給水溫度、給水壓力、蒸汽流量、蒸汽溫度、蒸汽壓力、二次側(cè)水位等。下降通道上的流量采用超聲波流量計進行測量。在下降段底部布置2個熱電偶用于測量進入管束區(qū)的流體溫度。

對于一次側(cè)沿程流體溫度，采用K型鎧裝熱電偶，熱電偶采用管接頭方式從傳熱管內(nèi)部引出，再采用管接頭將其引出筒體。安裝方式如圖5所示。管接頭作為鎧裝熱電偶穿過壓力邊界的密封元件，在與傳熱管焊接安裝前進行焊接工藝試驗，驗證了密封可靠性。

圖5 一次側(cè)測溫熱電偶管接頭安裝方式

空泡份額的測量分別采用光纖探針法、伽馬射線法和壓差法獲得局部點、截面和體積含汽率。

光纖探針法是利用光學的全反射原理，傳輸光在光纖探針敏感頭界面與媒質(zhì)界面處符合全反射條件時，這一部分光將返回敏感頭內(nèi)部，如果敏感頭接觸的媒質(zhì)折射率產(chǎn)生變化時，由于可能將不滿足全反射條件，此時不能產(chǎn)生全反射現(xiàn)象。通過檢測經(jīng)敏感頭反射回接收端的光強的強弱變化而產(chǎn)生的高低不同的電信號來分辨光纖探針敏感頭測量到的是氣相還是液相，從而根據(jù)產(chǎn)生的連續(xù)變化的電壓信號測量局部截面持氣率，光纖探針敏感頭測量原理圖如圖6所示。試驗中測點位置見圖7。具體測點位置和探針數(shù)量見表3。

射線法是一種非接觸式的測量方法，其原理是根據(jù)汽、液兩種物質(zhì)對射線的吸收程度不同推算出射線照射區(qū)域的含汽率。該方法不會破壞管道中的流場和溫度場的自然分布，適用于高溫高壓條件下測量，原理如圖8所示。試驗測點位置如圖9所示。

圖6 光纖探針法測量含汽率的原理示意

圖7 空泡份額測點布置

表3 光纖探針空泡份額測點位置

圖8 射線法原理示意

圖9 射線法測點位置示意

壓差法則是通過測量兩個取壓截面的壓力差，計算兩個截面之間的流體體積含汽率。

為了獲得較多的試驗數(shù)據(jù)，試驗工況圍繞額度工況，分別改變一次側(cè)進出口溫度、給水溫度、加熱功率和循環(huán)倍率等，完成了11個試驗工況。

試驗測量的一、二次側(cè)溫度分布如圖10所示。三種方法空泡份額測量結(jié)果見表4。

2.2 試驗件模擬分析

采用本文開發(fā)的程序?qū)υ囼灱M行數(shù)值模擬分析，預測二次流場數(shù)據(jù)。計算中選取試驗件管板以上至分離器入口為分析對象。對計算區(qū)域進行建模，管束區(qū)采用多孔介質(zhì)模型進行模擬，為了方便設置支承板阻力系數(shù)，支承板阻力采用多孔階躍邊界描述，在邊界上設定支承板的阻力系數(shù)參數(shù)和厚度。為了保證數(shù)值計算穩(wěn)定性與準確性，采用一體化網(wǎng)格，對整個計算區(qū)域進行建模，然后劃分網(wǎng)格，所劃分網(wǎng)格全為六面體網(wǎng)格，網(wǎng)格總數(shù)為196 107，計算模型網(wǎng)格如圖11所示。

(a)一次側(cè)溫度分布

(b)二次側(cè)溫度分布

表4 空泡份額測量結(jié)果對比

根據(jù)試驗件和試驗工況，計算中入口為下部開口，出口為分離器入口，入口邊界設為速度入口，不同工況給出不同入口流速，出口邊界設為壓力出口。計算中采用COUPLE算法求解控制方程組。本分析采用UDF輸入自定義模塊，如附加源項、各附加阻力項，因此計算過程中，除求解控制方程，還需調(diào)用UDF，在User-Defined Profile定義計算區(qū)域孔隙率；在User-Defined Init中定義熱源分布、網(wǎng)格潤濕面積，以及其他需要在計算進行之前計算的變量；在User-Defined Adjust中計算求解阻力源項所需的各個變量，如各控制體的阻力系數(shù)。

計算結(jié)果如圖12所示。為了圖片更加清晰，圖片在高度方向縮短了15倍。

圖11 計算模型網(wǎng)格

圖12 工況1的二次側(cè)溫度云圖及氣相體積份額云圖

2.3 程序?qū)Ρ闰炞C

為了驗證程序分析的準確性，將相同工況下程序預測結(jié)果與試驗測量結(jié)果進行比較。二次側(cè)溫度進行了歸一化處理，用測量點的溫度除以該工況下的飽和溫度，熱側(cè)和冷側(cè)歸一化比較結(jié)果如圖13所示?？张莘蓊~比較如圖14所示。

(a)熱側(cè)

(b)冷側(cè)

(a)射線法測量的空泡份額

(b)探針法測量的二次側(cè)空泡份額

3 結(jié)論

本文針對壓水堆核電廠蒸汽發(fā)生器結(jié)構(gòu)特點，以CFD軟件中多孔介質(zhì)模型為基礎，采用UDF對CFD軟件進行了二次開發(fā)，完成了蒸汽發(fā)生器二次側(cè)三維兩相熱工水力分析程序的開發(fā)。并基于?；治觯捎眯∫?guī)模管束進行流動傳熱試驗，模擬蒸汽發(fā)生器一、二次側(cè)流場，獲得試驗體試驗工況下的一、二次側(cè)流場數(shù)據(jù)，采用該程序?qū)υ囼灱７治觥＝Y(jié)果表明：程序能夠?qū)崿F(xiàn)對蒸汽發(fā)生器試驗件的模擬，預測結(jié)果與試驗測量結(jié)果符合良好，證明了該程序計算結(jié)果的可靠性。

未來蒸汽發(fā)生器二次側(cè)三維兩相熱工水力分析程序還需開展進一步工作。經(jīng)后續(xù)不斷改進完善，程序?qū)崿F(xiàn)在工程上的應用，該程序可成為我國核電廠蒸汽發(fā)生器設計研發(fā)的重要工具。