房 鵬 楊永偉 趙澤龍

1（中國科學(xué)院近代物理研究所 蘭州 730000）

2（中國科學(xué)院大學(xué)核科學(xué)與技術(shù)學(xué)院 北京 100049）

快堆（Fast Reactor）具有可增殖核燃料、嬗變長壽命放射性廢物等方面的優(yōu)勢。快堆消耗的是豐度高、不易裂變的238U。研究證明，相比壓水堆的0.5%～1%的利用率，功率較大的快堆可以將天然鈾利用率達(dá)到60%～70%［1］，同時還能嬗變錒系核素，減少污染物質(zhì)的排放?？於咽强沙掷m(xù)發(fā)展中重要的堆型，快堆的安全問題是核能發(fā)展中的重要問題，因此對于快堆來說安全問題的研究分析有實際工程意義。

在核反應(yīng)堆安全事故分析中，構(gòu)建精確的回路模型，就要求解大量的物理現(xiàn)象的偏微分方程組，會花費(fèi)大量時間和費(fèi)用，而且對于問題的初步分析和工作人員的監(jiān)督性分析是不方便的。如果想得到初步性的分析結(jié)果，或是事故過程中各種物理現(xiàn)象的變化圖像，集總參數(shù)法往往能提供足夠精確的答案［2］?？於逊譃榛芈肥脚c池式，池式快堆相比回路式具有冷卻劑喪失率低、冷池?zé)釕T性大、緩解事故能力強(qiáng)、結(jié)構(gòu)安全等優(yōu)點，因此池式堆應(yīng)用更廣泛［3］。本文以池式快堆（Pool-type Fast Reactor）建立集總參數(shù)模型并進(jìn)行無保護(hù)失流事故的模擬與分析，不同于其他反應(yīng)堆集總參數(shù)模型，本文建立了完整的反應(yīng)堆一回路模型并計算回路中堆芯與中間換熱器（Intermediate Heat eXchanger，IHX）的物理過程，并將結(jié)果與實驗和其他機(jī)構(gòu)計算值進(jìn)行比對，驗證了集總參數(shù)法在池式鈉冷快堆無保護(hù)失流事故安全分析中的可行性。

1 中子學(xué)計算

該程序中堆芯功率的計算采用點堆模型，6組緩發(fā)中子的點堆動力學(xué)方程如下：

點堆中子動力學(xué)模型是計算反應(yīng)堆動態(tài)學(xué)中最常用的方法，但因為耦合的點堆動力學(xué)微分方程組存在很強(qiáng)的剛性，所以需要采用隱式求解，該程序采用分段多項式法［4］進(jìn)行求解，該方法對解決強(qiáng)剛性的快堆問題更為有效［5-6］。

2 反應(yīng)堆回路模型及計算方法

反應(yīng)堆回路的模型化是簡化計算、使用集總參數(shù)法完成反應(yīng)堆瞬態(tài)分析的前提與基礎(chǔ)。根據(jù)池式堆的布置以及集總參數(shù)法的需要，將池式堆各重要部分簡化，建成完整的反應(yīng)堆一回路系統(tǒng)，其模型化如圖1所示。

圖1 反應(yīng)堆一回路模型Fig.1 Model of the reactor primary loop

堆芯的燃料組件與流經(jīng)堆芯的冷卻劑的傳熱過程，簡化為三個節(jié)點的傳熱：燃料棒、包殼、冷卻劑。圖2為堆芯傳熱計算的模型。

圖2 堆芯基礎(chǔ)計算模型Fig.2 Basic computational model of the reactor core

堆芯建模采用單通道模型，堆芯冷卻劑入口溫度采用下腔室冷卻劑的出口溫度。

中間換熱器采用立式布置的管殼式換熱器［7］，如圖3所示。

換熱器內(nèi)考慮了三個換熱節(jié)點：主回路冷卻劑、換熱器管道、次回路冷卻劑，如圖4所示。

堆芯內(nèi)燃料棒與包殼中的氣隙部分換熱過程忽略，簡化為燃料棒直接與包殼換熱，氣隙在其中的作用簡化為燃料棒與包殼之間的等效換熱系數(shù)。假設(shè)燃料棒與包殼始終不接觸，等效換熱系數(shù)計算方法［8］為：

式中：λg為氣體導(dǎo)熱系數(shù)；δg為間隙厚度；σ為斯蒂芬-玻爾茲曼系數(shù)；Σu和Σc為燃料和包殼的表面發(fā)射率；Tus為燃料元件芯塊外表面溫度；Tci為包殼內(nèi)表面溫度。式（6）中后一項為輻射傳熱項，將其忽略，保留氣體導(dǎo)熱項。即：

式中：MF為燃料質(zhì)量；CF為燃料比熱；P（t）為當(dāng)前時刻的功率；s為二者之間的換熱面積；TˉF為燃料平均溫度；TˉClad為包殼平均溫度。

包殼與冷卻劑之間傳熱的能量方程為：

式中：MP為主回路冷卻劑質(zhì)量；CP為冷卻劑比熱；hc為冷卻劑與包殼間的傳熱系數(shù)；Sc為二者傳熱面積；TˉC為堆芯冷卻劑平均溫度；W為主回路冷卻劑流量；TCin為堆芯冷卻劑入口溫度。

同時作為燃料與冷卻劑之間的節(jié)點，包殼的能量方程為：

圖3 中間換熱器模型Fig.3 IHX model

圖4 換熱器計算模型Fig.4 Computational model of the IHX

對于換熱器內(nèi)的換熱過程，也采用集總參數(shù)法來計算，主回路側(cè)冷卻劑先與管道這一節(jié)點進(jìn)行換熱，再由管道與次回路側(cè)冷卻劑進(jìn)行換熱。

對于主回路冷卻劑與管道，能量方程如下：

式中：Me為換熱器內(nèi)主回路冷卻劑質(zhì)量；ht為主回路冷卻劑與管道的傳熱系數(shù)；St為主回路冷卻劑與管道的傳熱面積；TˉPIHX為換熱器內(nèi)主回路冷卻劑平均溫度；TˉT為換熱器管道平均溫度；TPin為主回路側(cè)冷卻劑入口溫度。

管道的傳熱能量方程為：

式中：MT為管道材料質(zhì)量；CT為管道材料比熱；TˉS為次回路冷卻劑平均溫度；he為管道與次回路冷卻劑間的傳熱系數(shù)；Se為管道與次回路冷卻劑間的換熱面積。

管道與次回路側(cè)冷卻劑的傳熱能量方程為：

式中：MS為次回路冷卻劑質(zhì)量；CS為次回路冷卻劑比熱；W′為次回路冷卻劑流量；TˉS為換熱器次回路側(cè)冷卻劑平均溫度；TSin為換熱器次回路側(cè)冷卻劑入口溫度。

對于進(jìn)出口溫度與平均溫度在集總參數(shù)算法中的關(guān)系，以堆芯冷卻劑為例：

對于反應(yīng)堆中的腔室和冷池模型，因為在計算中主要關(guān)注其出口溫度，并假設(shè)絕熱邊界條件，故選用修正的完美混合模型［7］：

式中：ε是有效混合因子；C是冷卻劑比熱容；M是腔體質(zhì)量；Q是質(zhì)量流量；Tin是腔室的入口溫度。

其中：

式中：Qe為穩(wěn)態(tài)流量。

在池式快堆中一般用液態(tài)鈉、液態(tài)鉛鉍作為冷卻劑。傳熱能力不僅受流量變化的影響，冷卻劑溫度的變化會引起物性的變化，繼而影響傳熱系數(shù)的計算，因此需要計算冷卻劑的物性參數(shù)。

液態(tài)金屬冷卻劑鈉［11］和鉛鉍［12］的物性參數(shù)變化如表1所示。

液態(tài)鈉［10］與液態(tài)鉛鉍［13］努賽爾數(shù)計算采用表2的方法。

表1 液態(tài)鈉和鉛鉍的熱物性Table 1 Properties of sodium and LBE

表2 液態(tài)鈉和鉛鉍的努賽爾數(shù)Table 2 Nusselt of sodium and LBE

本文的管束計算模型選用六角陣列，以堆芯為例，表2中的P為相鄰燃料棒中心距；D為燃料棒直徑；Pe為佩克萊數(shù)。

由表2中的方法得到冷卻劑的Nu數(shù)，再根據(jù)式（17）得出傳熱系數(shù)：

式中：h為傳熱系數(shù)；L為特征長度；K為冷卻劑熱導(dǎo)率。

Pe在熱量傳輸中計算方法為Pe=Re·Pr；Re為雷諾數(shù)；Pr為普朗特數(shù)，二者計算公式分別為：

式中：ρ為冷卻劑密度；υ為冷卻劑流速；d為特征長度；η為冷卻劑黏度；C為冷卻劑比熱；K為冷卻劑熱導(dǎo)率。

3 EBR-II基準(zhǔn)題驗證

本文選取IAEA（International Atomic Energy Agency）提供的美國阿貢國家實驗室（Argonne National Laboratory，ANL）鈉冷快堆EBR-II的基準(zhǔn)題［10，14］，根據(jù)其 SHRT-45R（Shutdown Heat Removal Tests）實驗的結(jié)果進(jìn)行驗證。根據(jù)基準(zhǔn)題中的圖表確定其各部分幾何的相關(guān)參數(shù)。

3.1 實驗描述

在SHRT-45R開始之前，EBR-II在全功率全流量下運(yùn)行了足夠長的時間，使系統(tǒng)達(dá)到穩(wěn)態(tài)。在測試開始之前，控制棒的驅(qū)動器被停用。SHRT-45R的瞬態(tài)是通過打開每個泵的電機(jī)發(fā)電機(jī)組的電源斷路器，從而切斷泵的電源來啟動的。結(jié)果導(dǎo)致的情況與核電站停電時的情況相似。隨著SHRT-45R測試的繼續(xù)進(jìn)行，反應(yīng)堆功率由于測試開始后的反應(yīng)性反饋而下降。

該基準(zhǔn)題為核熱耦合問題，邊界條件為主、次回路流量的變化以及次回路冷卻劑入口溫度的變化。中子學(xué)計算部分主要考慮4種反應(yīng)性反饋：燃料溫度反饋、冷卻劑溫度反饋、堆芯徑向膨脹、燃料軸向膨脹。其參數(shù)參考ANL提供的數(shù)值，其中多普勒系數(shù)為自行計算，反應(yīng)性反饋系數(shù)如表3所示。

點堆方程計算所需要的六組緩發(fā)中子份額與先驅(qū)核衰變常數(shù)數(shù)值［7］如表4所示。

表3 反應(yīng)性反饋系數(shù)Table 3 Reactivity feedback coefficients

表4 緩發(fā)中子數(shù)值Table 4 Delayed neutron data

3.2 穩(wěn)態(tài)計算

穩(wěn)態(tài)的計算是進(jìn)行瞬態(tài)分析的前提，依照基準(zhǔn)題［10］所提供的數(shù)據(jù)對穩(wěn)態(tài)進(jìn)行計算，與基準(zhǔn)題進(jìn)行對比。

表5為程序計算的穩(wěn)態(tài)數(shù)值與實驗值的對比，對比結(jié)果可知穩(wěn)態(tài)計算值與實驗值吻合較好，可以在此基礎(chǔ)上進(jìn)行瞬態(tài)計算。

表5 穩(wěn)態(tài)計算值Table 5 Steady-state calculation

3.3 瞬態(tài)分析

瞬態(tài)的邊界條件由基準(zhǔn)題［10］中提供的數(shù)據(jù)與圖像給出，包括了一回路流量、換熱器二回路入口溫度、二回路流量隨時間的變化。

圖5 回路中的流量變化Fig.5 Flow rate changes in primary loop and secondary loop

圖6 換熱器二次側(cè)入口溫度Fig.6 Inlet temperature at secondary side of IHX

圖7 功率計算結(jié)果對比Fig.7 Comparison of power calculation results and experimental results

這里強(qiáng)迫對流循環(huán)過渡到自然對流循環(huán)的過程中，一回路流量由基準(zhǔn)題［10］提供的主泵轉(zhuǎn)速數(shù)據(jù)按照轉(zhuǎn)速流量正比關(guān)系求得。二回路流量變化曲線由基準(zhǔn)題直接給出。

圖6數(shù)據(jù)由基準(zhǔn)題［10］直接給出，該數(shù)據(jù)為式（13）中換熱器次回路側(cè)冷卻劑入口溫度TSin的變化，作為換熱器瞬態(tài)計算的邊界條件。

圖7為程序計算的瞬態(tài)過程中功率變化與實驗結(jié)果的對比。

瞬態(tài)過程中前100 s功率下降略快，100 s后計算值偏高，考慮反饋系數(shù)的選取與平均溫度的計算與實驗存在差別造成的結(jié)果誤差。對比其他研究機(jī)構(gòu)結(jié)果［14］該結(jié)果在可接受范圍內(nèi)。

圖8和圖9分別為堆芯內(nèi)燃料平均溫度和冷卻劑平均溫度變化以及4種主要反應(yīng)性反饋的變化情況。

圖8 燃料與冷卻劑平均溫度Fig.8 Average temperatures of fuel and coolant

圖9 主要反應(yīng)性反饋計算結(jié)果Fig.9 Results of feedback computation of primary reactivity

圖10為堆芯冷卻劑入口溫度變化實驗值與程序計算值的對比。

圖10 堆芯入口溫度變化對比Fig.10 Comparison of temperature changes at core inlet

堆芯入口溫度變化與實驗結(jié)果有所差別是因為鈉池的計算采用了完美混合模型。實際上因為冷池的存在，堆芯入口溫度在整個瞬態(tài)過程中變化不明顯，且計算值對比實驗值誤差較小，對后續(xù)計算不會產(chǎn)生較大影響，所以認(rèn)為堆芯入口溫度的計算是可靠的。

Z型管是一條連接堆芯上腔室和換熱器入口的管道，Z型管入口溫度表征了堆芯出口溫度，因此是一項重要考察指標(biāo)。圖11為程序計算的Z型管入口溫度變化情況與實驗值和其他機(jī)構(gòu)計算值［14］的比較。因為實驗數(shù)據(jù)的部分缺失，所以導(dǎo)致實驗值的變化曲線存在間斷。圖11中機(jī)構(gòu)對上腔室處理方法有多種：西安交通大學(xué)（Xi’an Jiaotong University，XJTU）與福井大學(xué)（University of FUKUI）采用管道模型模擬上腔室，荷蘭核研究與咨詢集團(tuán)（Nuclear Research and Consultancy Group，NRG）采用了多控制體模型，韓國原子能研究院（Korea Atomic Energy Research Institute，KAERI）采用了熱結(jié)構(gòu)模型，而阿貢國家實驗室（Argonne National Laboratory，ANL）采用的是零維體積建模法。對比結(jié)果顯示：本文計算的Z型管入口溫度整體變化趨勢為先上升后下降，與實驗值變化一致，初期計算值上升速度較快，是因為本程序的上腔室采用了完美混合模型，沒有考慮上腔室中較熱流體流動引起的熱延時。末態(tài)溫度計算值偏高一方面是因為功率計算值偏高，使得堆芯冷卻劑溫度偏高；另一方面因為假設(shè)了絕熱邊界條件，沒有考慮上腔室與冷池之間的換熱。75～200 s內(nèi)實驗數(shù)據(jù)缺失，無法與實驗值進(jìn)行對比，因此峰值溫度選取與實驗結(jié)果最為接近的ANL計算值進(jìn)行比對，幾家機(jī)構(gòu)計算結(jié)果顯示，峰值溫度相對誤差最大為5.15%，最小為0.86%；末態(tài)溫度相對誤差最大為2.93%，最小為0.98%；本文峰值溫度相對誤差為4.05%，末態(tài)溫度相對誤差為2.37%，與其他機(jī)構(gòu)處于同一范圍內(nèi)，因此認(rèn)為程序計算結(jié)果可靠。

圖11 Z型管入口溫度結(jié)果對比Fig.11 Comparison of calculation results of Z-pipe inlet temperatures

4 結(jié)語

本文中，對池式快堆EBR-II建模并采用集總參數(shù)計算方法對SHRT-45R無保護(hù)失流工況進(jìn)行了模擬。計算結(jié)果表明：采用點堆模型與集總參數(shù)模型模擬得到的結(jié)果絕大部分與實驗結(jié)果較為相符，具有較高的準(zhǔn)確性，對比國內(nèi)外程序的計算結(jié)果也處在同一水平，同時大大節(jié)省計算時間成本與程序編寫工程量。因此，點堆模型結(jié)合集總參數(shù)的安全分析方法適用于池式鈉冷快堆無保護(hù)失流工況的模擬，并且該程序可對池式鈉冷快堆在事故中的堆芯行為與固有安全性做出較為準(zhǔn)確的預(yù)估計算。