明 楊，易經(jīng)緯，方華偉，劉 凱，趙富龍,*，譚思超，田瑞峰

(1.哈爾濱工程大學 核安全與仿真技術國防重點學科實驗室，黑龍江 哈爾濱 150001;2.中國核動力研究設計院，四川 成都 610213)

布雷頓循環(huán)反應堆是第4代核能系統(tǒng)的優(yōu)先發(fā)展堆型之一[1]。該堆型通常采用化學性質(zhì)穩(wěn)定、熱工性能良好的氣體冷卻劑(如氦氣、氙氣等)，堆芯冷卻劑出口溫度可達1 000 ℃甚至更高[2]，其主要特點是固有安全性高、經(jīng)濟性好、發(fā)電效率高、工藝熱應用廣泛[3]。

布雷頓循環(huán)系統(tǒng)分析程序方面，國內(nèi)外已開展了很多研究。成利等[4]利用simulink程序，對核電站反應堆溫度控制系統(tǒng)進行建模，驗證動態(tài)系統(tǒng)的建模方法可行性。李海鵬等[5]采用集總參數(shù)的建模方法，建立了HTR-10系統(tǒng)的動態(tài)模型，結果能正確反映系統(tǒng)的正常運行特性。唐麗麗[6]針對燃氣輪機布雷頓循環(huán)進行了建模仿真，研究了壓氣機與氣輪機的運行特性。王恩華等[7]針對布雷頓循環(huán)氦氣輪機的起動過程進行了動態(tài)仿真，得到了壓氣機和透平的工作曲線。解衡等[8]開發(fā)了仿真程序，對布雷頓循環(huán)負荷喪失事故工況以及主換熱器換熱量減少兩種瞬態(tài)工況進行了仿真。李文龍等[9]開發(fā)出布雷頓循環(huán)的瞬態(tài)分析軟件，針對引入反應性、工質(zhì)質(zhì)量流量變化等工況進行了仿真研究。Wright等[10-11]基于simulink程序首先針對閉式布雷頓循環(huán)建立了系統(tǒng)分析程序，分析系統(tǒng)的變工況特性，并在已有程序基礎上引入控制，研究不同控制方法下系統(tǒng)的響應特性。Verkerk[12]針對布雷頓循環(huán)，預測了布雷頓循環(huán)在反應堆變功率工況下的動態(tài)變化。Becker等[13]采用計算流體動力學軟件CFX-4建立了高溫氣冷堆傳熱流動的三維數(shù)值仿真模型，進行了穩(wěn)態(tài)和瞬態(tài)仿真計算與分析。Tauveron等[14-15]開發(fā)出一維計算機仿真程序CATHARE對布雷頓循環(huán)進行安全性分析，并對系統(tǒng)破口事故進行了仿真。Saez等[16]利用CATHARE程序?qū)Σ祭最D循環(huán)的氣輪機組通流阻塞情況進行了仿真研究。Bentivoglio等[17]應用CATHARE程序?qū)Φ聡鳲berhausenⅡ50 MW直接布雷頓循環(huán)核電站進行了仿真研究。

調(diào)研發(fā)現(xiàn)，國內(nèi)外學者進行了諸多關于布雷頓循環(huán)反應堆系統(tǒng)的仿真分析研究，并進行了系統(tǒng)運行特性的分析以及參數(shù)優(yōu)化，但針對小型直接布雷頓循環(huán)反應堆系統(tǒng)分析程序的相關研究較少。此外工質(zhì)物性模型多為常物性，對系統(tǒng)變工況運行時的仿真能力存在一定的不足[10-11]。因此本文針對小型直接布雷頓循環(huán)反應堆系統(tǒng)建立系統(tǒng)分析程序，并建立精確的工質(zhì)物性數(shù)據(jù)庫，考慮工質(zhì)變工況時的物性變化，開展小型布雷頓循環(huán)系統(tǒng)穩(wěn)態(tài)工況仿真以及反應堆升降功率、引入階躍反應性、關鍵設備參數(shù)變化、機組甩負荷等瞬態(tài)工況下的動態(tài)仿真。

1 系統(tǒng)分析程序建模

1.1 建模方法及假設

建模對象為一種以氦氙混合氣體作為工質(zhì)的小型氣冷快堆，系統(tǒng)循環(huán)方式為直接布雷頓循環(huán)，其原理如圖1所示，相應的一維系統(tǒng)分析程序基于圖1建立，子模塊功能和布置與系統(tǒng)圖基本相同，旨在盡可能還原系統(tǒng)的真實運行過程。

圖1 直接布雷頓循環(huán)氣冷反應堆系統(tǒng)原理圖

基本建模思路為：1) 依據(jù)布雷頓循環(huán)各重要設備的熱力學特性、流體力學特性、轉(zhuǎn)動特性等逐步建立各設備的數(shù)學物理模型，基于質(zhì)量守恒、能量守恒、轉(zhuǎn)子功率平衡等原理，構造各重要設備共同運行的數(shù)學物理關系；2) 通過數(shù)值方法對模型進行求解，獲得主要設備不同工況下的運行數(shù)據(jù)，模擬布雷頓循環(huán)的實際工作過程，并加入相應的控制方法，盡可能復現(xiàn)布雷頓循環(huán)的實際工作狀態(tài)。

建立的系統(tǒng)分析程序包括8個子模塊：反應堆模塊、氣輪機模塊、壓氣機模塊、回熱器模塊、冷卻器模塊、轉(zhuǎn)子模塊、流量變化模塊、PID控制模塊。

各子模塊內(nèi)部的計算公式、自定義函數(shù)、物性數(shù)據(jù)庫、曲線插值方法等實現(xiàn)了每個子模塊的計算功能，選取合適的初始參數(shù)與求解算法，進而完成布雷頓循環(huán)的系統(tǒng)分析計算。為簡化數(shù)學物理模型，基本假設如下：1) 工質(zhì)為實際氣體，在關鍵設備內(nèi)一維流動，忽略工質(zhì)物性在空間上的不均勻性；2) 忽略工質(zhì)泄漏，系統(tǒng)在運行時與外界環(huán)境無質(zhì)量與能量的交換；3) 忽略轉(zhuǎn)動部件運行時摩擦導致的能量損失；4)忽略長期運行對換熱系數(shù)、局部阻力系數(shù)等參數(shù)的影響。

1.2 關鍵設備物理模型

1) 反應堆模型

反應堆是布雷頓循環(huán)的能量來源，反應堆模塊需在各工況下正確計算反應堆的進出口工質(zhì)溫度、壓力以及反應堆功率，同時反映反應堆的溫度反饋與核反饋等物理過程。由于在系統(tǒng)動態(tài)過程中，主要關注反應堆功率變化以及中子通量變化，因此選取點堆中子動力學方程，主要方程如下：

(1)

(2)

(3)

ρ=ρin-α(Tf-Tf0)

(4)

式中：N堆為反應堆熱功率；ρ為總反應性；β為總緩發(fā)中子份額；βi為第i組緩發(fā)中子份額；Λ為中子代時間；λi為第i代緩發(fā)中子衰變常量；Ci為緩發(fā)中子先驅(qū)核濃度；ρin為控制棒引入的反應性；α為溫度負反饋系數(shù)；Tf為核燃料溫度；Tf0為額定功率下燃料溫度。

堆芯傳熱模型基于集總參數(shù)法建立，由于堆芯傳熱模型關心反應堆的出口溫度與堆芯功率的關系，可認為整個堆芯的燃料元件在同一時刻均處于同一溫度。模型主要方程如下：

(5)

G堆cp冷(T堆冷o-T堆冷i)]

(6)

(7)

式中：Pr為反應堆熱功率；Tf為燃料平均溫度；ρf為燃料元件密度；cvf為燃料元件熱容；Vc為反應堆總體積；α為燃料元件換熱系數(shù)；A為燃料元件換熱面積；ρ冷為冷卻劑密度；cp冷為冷卻劑比定壓熱容；T堆冷a為堆芯冷卻劑平均溫度；T堆冷o為堆芯冷卻劑出口溫度；T堆冷i為堆芯冷卻劑進口溫度；ε為堆芯孔隙率；G堆為堆芯質(zhì)量流量。

2) 氣輪機與壓氣機模型

由于系統(tǒng)需計算氣輪機與壓氣機的變工況運行特性，準穩(wěn)態(tài)法不適用于偏離穩(wěn)態(tài)工況較大的計算，因此系統(tǒng)分析程序采取特性曲線插值法對氣輪機進行建模。氣輪機模型主要方程如下：

δ氣=T氣i/T氣o

(8)

(9)

(10)

式中：G氣為氣輪機質(zhì)量流量；T氣o為氣輪機出口溫度；T氣i為氣輪機進口溫度；δ氣為氣輪機膨脹比；n為氣輪機轉(zhuǎn)速；η氣為氣輪機等熵效率；N氣為氣輪機功率；k為絕熱指數(shù)。

系統(tǒng)分析程序采取特性曲線插值法對壓氣機進行建模。壓氣機模型主要方程如下：

ε壓=T壓o/T壓i

(11)

(12)

(13)

式中：G壓為壓氣機質(zhì)量流量；T壓o為壓氣機出口溫度；T壓i為壓氣機進口溫度；ε壓為壓氣機壓縮比；η壓為壓氣機等熵效率；N壓為壓氣機功率；k為絕熱指數(shù)。

3) 回熱器與冷卻器模型

回熱器模型基于其結構和換熱方式對其建模，并作一定程度的假設，忽略回熱器與外界的換熱，認為同一截面內(nèi)工質(zhì)的溫度、速度、壓力參數(shù)一致，流體沿軸向一維流動，主要方程如下:

G低i=G低o=G壓

(14)

G高i=G高o=G氣

(15)

(16)

T高s=(T高o-T高i)/2

(17)

(18)

(19)

式中：G高i為回熱器高溫側進口質(zhì)量流量；G高o為回熱器高溫側出口質(zhì)量流量；T高i為回熱器高溫側進口溫度；T高o為回熱器高溫側出口溫度；G低i為回熱器低溫側進口質(zhì)量流量；G低o為回熱器低溫側出口質(zhì)量流量；T低i為回熱器低溫側進口溫度；T低o為回熱器低溫側出口溫度；Q高為回熱器高溫側換熱量；Q低為回熱器低溫側換熱量；M高為回熱器高溫側流體總質(zhì)量；M低為回熱器低溫側流體總質(zhì)量；a高、a低為回熱器換熱系數(shù)；A高、A低為回熱器換熱面積；e高、e低為回熱器換熱面效率。冷卻器模型與回熱器模型類似。

4) 轉(zhuǎn)子模型

轉(zhuǎn)子模型建立基于假設3，忽略轉(zhuǎn)動部件間的摩擦。系統(tǒng)處于穩(wěn)態(tài)運行時，轉(zhuǎn)子力矩平衡，氣輪機轉(zhuǎn)速維持恒定。當氣輪機的輸出功與壓氣機耗功不等時，轉(zhuǎn)子將在力矩的作用下加速或減速。轉(zhuǎn)子模型主要方程如下：

(20)

n=30ω/π

(21)

式中：ω為轉(zhuǎn)子角速度；Ne為系統(tǒng)電功率；Nt、Nc為氣輪機、壓氣機功率；I為轉(zhuǎn)子轉(zhuǎn)動慣量；n為壓氣機轉(zhuǎn)速。

5) 壓降模型

氣體冷卻劑在管道中的壓降采用沿程摩擦壓降進行計算，文中采取經(jīng)典的Darcy公式，選取合理的摩擦系數(shù)，確定流道的水力直徑、通道長度等參數(shù)即可計算管道中的壓降與關鍵設備中的壓力變化情況，主要方程如下：

(22)

式中：Δpf為壓降；f為摩擦系數(shù)，它與流體的流動性質(zhì)、流動狀態(tài)、受熱情況、流道的幾何形狀、表面粗糙程度相關；ρ為工質(zhì)密度；L為通道長度；De為流道水力直徑；V為體積。

6) 物性模型

首先查詢工質(zhì)在系統(tǒng)運行溫度、壓力范圍內(nèi)的絕熱系數(shù)、比熱容、密度，并在系統(tǒng)分析程序中生成數(shù)據(jù)庫。之后定義插值函數(shù)的計算方法，在計算過程的每個時間步長內(nèi)實時更新物性參數(shù)。

2 計算結果分析

2.1 穩(wěn)態(tài)計算結果分析

對布雷頓循環(huán)系統(tǒng)在額定工況下穩(wěn)態(tài)運行進行計算，并將計算結果與系統(tǒng)設計值進行對比后求出仿真誤差，結果列于表1。可看出在布雷頓循環(huán)系統(tǒng)額定工況的仿真過程中最大仿真誤差為6.38%，計算結果較為合理，說明開發(fā)的系統(tǒng)分析程序可靠，可繼續(xù)進行布雷頓循環(huán)系統(tǒng)變工況條件下的動態(tài)分析。

2.2 瞬態(tài)計算結果分析

1) 引入反應性工況

在布雷頓循環(huán)的變工況運行過程中，反應堆會由于控制棒的動作或其他原因而引入反應性，會導致反應堆的功率發(fā)生變化，從而引發(fā)整個系統(tǒng)運行參數(shù)的較大波動，因此引入反應性是探究系統(tǒng)特性的重要工況之一。計算初始反應堆在設計工況下運行，1 000 s時引入100 pcm的正反應性，而后在2 000 s時引入100 pcm的負反應性，其他輸入變量保持額定運行值不變，對系統(tǒng)各關鍵參數(shù)進行仿真，總仿真時長為3 000 s。

反應堆出口溫度、核燃料溫度的動態(tài)變化如圖2所示。引入正反應性導致中子通量與反應堆熱功率迅速升高，導致堆芯核燃料與反射層溫度上升。由于反應堆的核反饋與溫度反饋作用，等同于向反應堆引入了負反應性，以此來抑制核燃料和反射層溫度的上升，工質(zhì)出口溫度與核燃料的響應曲線會呈波動變化并逐漸平穩(wěn)。引入負反應性時與引入正反應性時功率趨勢相反。

圖3為關鍵設備功率曲線。由圖3a可看出，引入正反應性時，反應堆功率峰值達到11.8 MW，隨后在燃料溫度效應、慢化劑溫度效應等物理現(xiàn)象的影響下出現(xiàn)負反饋，并迅速下降。之后呈現(xiàn)功率高低震蕩現(xiàn)象，并最終趨于穩(wěn)定，功率提升至10.5 MW。引入負反應性時，反應堆功率谷值達到8.0 MW，隨后在反應堆的自穩(wěn)調(diào)節(jié)作用下出現(xiàn)回升，之后呈現(xiàn)功率高低震蕩現(xiàn)象后穩(wěn)定，功率降至9.6 MW。由圖3b、c、d可看出，引入正反應性后反應堆功率的上升導致工質(zhì)出口溫度升高，同時氣輪機做功功率提升，從而使轉(zhuǎn)子轉(zhuǎn)速提升。同時壓氣機功率也相應升高，發(fā)電機功率由額定工況下的2.63 MW增至2.71 MW，布雷頓循環(huán)系統(tǒng)的功率輸出能力提升。

表1 額定工況仿真結果

圖2 反應堆出口溫度與核燃料溫度曲線

圖3 關鍵設備功率曲線

2) 溫度階躍變化工況

初始系統(tǒng)運行穩(wěn)定，在2 000 s時堆芯進口工質(zhì)溫度由604.9 K階躍上升至620 K，其他輸入變量保持額定值不變，對系統(tǒng)進行仿真，總仿真時長為4 000 s。反應堆出口溫度、反應堆核燃料溫度和反應堆功率的計算結果如圖4所示。

從圖4可知，反應堆進口工質(zhì)溫度階躍升高，由于工質(zhì)與堆芯以及反射層間的換熱量下降，堆芯燃料溫度與反射層溫度上升，由于堆芯溫度具有負溫度反應性系數(shù)，因此引入負反應性，導致反應堆功率下降，進而使反應堆功率降低，同時堆芯燃料溫度也進一步降低。反應堆功率、堆芯燃料溫度和反應堆出口溫度的響應曲線會呈先升高后降低的波動變化，并逐漸趨于穩(wěn)定。

3) 流量階躍變化工況

在2 000 s時反應堆進口工質(zhì)質(zhì)量流量由 4.60 kg/s階躍降至4.50 kg/s，其他輸入變量保持額定值不變，對系統(tǒng)進行仿真，總仿真時長為4 000 s。反應堆出口溫度、核燃料溫度和反應堆功率仿真結果如圖5所示。反應堆進口氦氣質(zhì)量流量階躍下降，工質(zhì)與堆芯及反射層間的換熱量下降，導致堆芯燃料溫度與堆芯反射層溫度上升。由于堆芯核燃料具有負溫度反應性系數(shù)，溫度升高導致引入負反應性，導致反應堆功率下降，進而導致反應堆出口溫度與核燃料溫度的回落。反應堆出口溫度、核燃料溫度和反應堆功率的響應曲線會呈先上升后下降的波動變化，并逐漸趨于穩(wěn)定。

圖4 溫度階躍變化工況下反應堆出口溫度、核燃料溫度、反應堆功率曲線

圖5 流量階躍變化工況下反應堆出口溫度、核燃料溫度、反應堆功率曲線

4) 氣輪機組甩負荷工況

在布雷頓循環(huán)的變工況運行過程中，可能出現(xiàn)氣輪機組甩負荷情況，即發(fā)電機轉(zhuǎn)子與氣輪機組轉(zhuǎn)子脫扣，氣輪機組失去載荷的工況，系統(tǒng)運行參數(shù)會因此產(chǎn)生較大的波動。因此本節(jié)將著重研究甩負荷工況下的系統(tǒng)響應特性。初始系統(tǒng)運行穩(wěn)定，在2 000 s時堆芯進口工質(zhì)溫度由604.9 K階躍升至 620 K，其他輸入變量保持額定值不變，對系統(tǒng)進行仿真，總仿真時長為4 000 s。

轉(zhuǎn)子轉(zhuǎn)速、反應堆功率和反應堆出口溫度的計算結果如圖6所示。由于氣輪機組失去載荷，與發(fā)電機脫扣，導致轉(zhuǎn)子轉(zhuǎn)速由15 000 r/min迅速飆升至16 200 r/min。系統(tǒng)工質(zhì)質(zhì)量流量迅速增大，導致冷卻劑從反應堆帶走的熱量增加。反應堆出口溫度在流量的增加作用下迅速降低，同時由于溫度反饋的作用，反應堆功率增高。這種反應堆功率的迅速升高對布雷頓循環(huán)系統(tǒng)的影響十分嚴重，可能導致嚴重事故的發(fā)生。由程序計算結果顯示，系統(tǒng)對轉(zhuǎn)子引入轉(zhuǎn)速控制是極其必要的。

加入轉(zhuǎn)速控制后的轉(zhuǎn)子轉(zhuǎn)速、反應堆功率和反應堆出口溫度的計算結果如圖7所示?？刂葡到y(tǒng)相當于為氣輪機組加入引入的負載，當轉(zhuǎn)子轉(zhuǎn)速超過限值后對轉(zhuǎn)子進行剎車從而減少轉(zhuǎn)速，將轉(zhuǎn)子減速到15 370 r/min，從而抑制系統(tǒng)質(zhì)量流量的增加，并且控制反應堆功率在11 MW左右，從而控制反應堆功率在安全范圍內(nèi)，防止極限工況與事故發(fā)生。

圖6 轉(zhuǎn)子轉(zhuǎn)速、反應堆功率、反應堆出口溫度曲線

圖7 加入轉(zhuǎn)速控制后轉(zhuǎn)子轉(zhuǎn)速、反應堆功率、反應堆出口溫度曲線

3 結論

本文基于模塊化建模思想，建立了直接布雷頓循環(huán)反應堆系統(tǒng)分析程序，并針對典型穩(wěn)態(tài)與瞬態(tài)工況條件下系統(tǒng)的動態(tài)特性進行了模擬仿真，主要得到以下結論：

1) 系統(tǒng)分析程序穩(wěn)態(tài)仿真精度較高，關鍵設備參數(shù)仿真值與設計值相對誤差不超過6.38%；

2) 各子模塊耦合性能較好，能實時反映變工況下系統(tǒng)主要設備相關參數(shù)的變化趨勢；

3) 系統(tǒng)具有一定的自穩(wěn)自調(diào)能力，與小型布雷頓循環(huán)反應堆系統(tǒng)的固有安全性吻合較好；

4) 針對氣輪機組甩負荷這種極限工況，在引入控制模塊后，系統(tǒng)能在安全范圍內(nèi)運行，并驗證了引入控制的必要性。

綜上所述，本文開發(fā)的小型系統(tǒng)分析程序模塊化程度高、計算速度快，為系統(tǒng)運行提供分析與變工況特性預測，并為直接布雷頓循環(huán)氣冷反應堆系統(tǒng)的分析和優(yōu)化提供技術支撐。