快中子通量實驗堆失流事故三維數值模擬

周 磊，張大林，劉雅鵬，王式保，王心安，王成龍，田文喜，秋穗正，蘇光輝

(西安交通大學 動力工程及多相流國家重點實驗室，核科學與技術學院，陜西 西安 710049)

為準確評估鈉冷快堆在各類事故工況下的瞬態響應以及各類參數能否滿足事故工況下的驗收準則，通常采用系統程序針對鈉冷快堆進行建模計算以獲取事故工況下鈉冷快堆的各種關鍵參數，如西安交通大學開發的THACS[1-2]、中國原子能科學研究院開發的FR-Sdaso[3]和FRSYS[4]、華北電力大學開發的SAC-CFR[5-6]，以及美國阿貢國家實驗室開發的SAS4A/SASSYS-1[7]、俄羅斯國家科學中心物理動力研究院開發的GRIF[8]、法國原子能委員會開發的CATHARE[9]等，這類系統程序通常采用集總參數法或一維分析方法，以及有限的二維和偽三維分析方法以分析各類工況下的各種現象，為保證計算結果的可靠性，通常需要大量的實驗數據用以提供各種經驗關系式或物理模型。該類系統程序雖能滿足鈉冷快堆的設計和安全分析，但開發驗證成本高，適用范圍相對窄，同時在一定程度上忽略了反應堆系統局部現象，犧牲了核動力系統的經濟性。伴隨著高性能計算機的發展，對于核動力系統中部分系統與部件采用的全三維模擬逐漸成為可能。相較于系統程序，三維CFD方法能有效捕捉反應堆復雜的三維流動換熱現象，識別潛在的安全影響因素，充分挖掘反應堆的經濟性[10]。但不同于系統程序針對特定反應堆長期開展的大量驗證與確認工作以及采用偏保守的物理模型，三維CFD方法從物理模型上、數值算法上具有更加廣泛的適用范圍與適用對象，與真實的物理現象更為接近，但針對不同流態(如層流、湍流、過渡流)、不同工質(如水、液態鈉、液態鉛、熔鹽等)、不同工況(如強迫循環、自然循環工況)、不同流動現象(如回流、射流、繞流等)、不同流型(如單相流、氣液兩相流、氣固兩相流等)，CFD方法仍缺乏普適的湍流模型和數值算法，特別是在針對反應堆全堆級別的事故工況下的瞬態模擬驗證結果仍有限。目前已有的全堆數值模擬[11-12]采用商業軟件針對實驗快堆的失給水工況和主泵卡死等事故工況進行全三維模擬，但建模時對堆芯組件結構進行了較大的簡化，無法評估鈉冷快堆重要的盒間流現象[13]在瞬態過程中對反應堆安全性能的影響，此外模擬結果缺乏足夠的實驗數據的驗證，另一種較通用的方法是采用跨維度耦合[14-15]的方式進行模擬，針對核心關注部分采用三維CFD方法進行模擬，而針對其余部分采用系統程序進行模擬，該種方法在一定程度上結合了系統程序與CFD方法的優勢，保證了精度，但對于兩種不同程序的耦合又會帶來收斂問題，同時對于瞬態變化劇烈的工況，難以獲取較好的模擬結果。

因此，為確保CFD方法在針對鈉冷快堆的低Pr流體液態金屬鈉[16]，以及鈉冷快堆自然循環工況下[17]計算的準確性，本文利用快中子通量實驗堆(FFTF)進行失流事故實驗，并針對其建立相應的三維模型，其中冷池和熱池按照原型尺寸建模，堆芯組件簡化為多孔介質，同時保留堆芯盒間幾何特性，針對無保護失流事故進行900 s的瞬態模擬，并將模擬結果與實驗數據進行對比，初步驗證CFD方法在針對鈉冷快堆自然循環工況模擬中的可行性與可靠性。

1 FFTF系統及失流事故簡介

FFTF[18]是由美國西屋公司設計建造的一座400 MW熱功率，采用氧化物燃料、液態金屬鈉冷卻的實驗堆。FFTF位于華盛頓漢福德研究中心，于1980年首次臨界，一直運行到1992年，并于1993年關閉，期間進行了大量實驗與測試，為鈉冷快堆的燃料、材料、部件積累了大量經驗。1986年，為測試FFTF非能動安全特性，在該堆上進行了大量無保護瞬態實驗，其中包含了13次無緊急停堆的失流事故實驗(LOFWOS)，為液態金屬堆的程序驗證、鈉冷快堆的固有安全性提供了重要支撐。2018年，在IAEA的組織下，選取LOFWOS#13作為國際基準題，多個國家和組織參與了此次基準題計算，西安交通大學作為參與單位之一，開展了系統程序和三維瞬態計算。

圖1 冷卻系統示意圖Fig.1 Coolant system overview

FFTF是一座環路式鈉冷快堆，其典型環路及主要部件示于圖1。FFTF包含2個回路、3個環路，液態鈉從反應堆容器離開，分別進入3個熱腿，在經過每條環路上的中間熱交換器時將熱量傳遞給二回路中無放射性的鈉，最后通過二回路的空氣熱交換器將熱量導入環境中，因此該堆并不直接用來發電。FFTF還采用了獨特的安全設計以保證反應堆的安全，如較為特殊的氣體膨脹組件，能減小無保護失流事故的事故后果，該組件插入到堆芯徑向反射區的靠近內堆芯區，如圖2所示，組件內部頂端通過不銹鋼插件限制了約0.028 3 m3的氬氣，造成頂端氬氣無法泄漏，底端與柵板聯箱連通，液鈉可進入組件。因此，組件的液鈉高度與入口鈉的壓力及氬氣溫度有關，額定工況下，由于主泵提供較大的揚程，使得組件內液鈉的高度高于堆芯活性區，而當堆芯溫度升高或失去強迫循環后，液鈉的高度下降，堆芯徑向中子泄漏增加，為反應堆提供負反饋，從而保障反應堆的安全性。堆芯組件分布示于圖3，圖中DF為燃料組件，CR、SR為控制棒組件，GEM為氣體膨脹組件。圖3中具有特殊標記的組件PIOTA為裝有溫度測點的實驗組件，LOFWOS#13實驗開始前保持50%額定功率、100%額定流量、約7.5 h，使得整個系統達到熱平衡狀態。為保證實驗完整運行，反應堆安全系統在實驗前進行了調整，避免控制棒過早插入，對堆內流量和通量反應堆關閉信號也進行了關閉或旁通。實驗開始的初始時刻，主回路的3個環路主泵同時關閉，轉入惰轉，約90 s后主泵惰轉結束，反應堆轉入自然循環，整個實驗過程中二回路泵仍保持額定運行狀態。在實驗進行過程中未進行其他動作與干預。實驗開始前反應堆穩態參數列于表1。

圖2 氣體膨脹組件Fig.2 Gas expansion module

圖3 堆芯組件類型及分布Fig.3 Assembly type and distribution in core

表1 反應堆穩態參數Table 1 Steady state parameter of reactor

1.1 網格模型和計算邊界

FFTF堆本體如圖4所示，堆內構件眾多，結構復雜，尺寸跨度大，完全按照原始幾何建模劃分網格，將導致計算量巨大，現階段難以完成瞬態計算。因此，針對真實反應堆結構進行了一定程度的簡化。本文建模中忽略了一系列與主要流道不相關的幾何腔室與結構，如圖5所示，僅保留了與流動密切相關的幾何結構，如冷池、熱池、堆芯等關鍵結構，此外將各分隔板簡化為無厚度的面，對于位于柵板聯箱中的堆芯支持結構以及入口管腳，該結構復雜且數目眾多，主要起流量分配的作用，因此在幾何建模過程中忽略了該部分結構，并將該部分結構的壓降及流量分配的作用歸到燃料組件的入口處。堆芯的各類組件采用繞絲棒束組件，流道狹長，針對每盒組件進行了多孔介質簡化，最終的網格模型如圖6所示，其中堆芯組件采用六面體網格，其余部分采用四面體網格，整體網格呈混合網格，各類型網格面連續一致。

圖4 FFTF堆本體Fig.4 Overview of reactor vessel

圖5 FFTF簡化幾何模型Fig.5 Simplified geometry model of FFTF

1.2 數學模型與邊界條件

1) 堆芯組件模型

由于對反應堆堆內結構進行了一定程度的簡化，特別是堆芯中忽略了組件插入柵板聯箱中的組件分流結構，可能直接導致堆芯各類組件的流量分配與設計值產生較大差距；另外堆芯各類幾何結構復雜，形狀阻力難以確定，為保證簡化后流量分配與設計值保持一致，在穩態下采用如下方法進行組件進口區域的多孔介質參數調整。

圖6 堆本體網格模型Fig.6 Mesh model for reactor vessel

將堆芯內各組件按照各組件的實際結構沿軸向標記為進口區、出口區、繞絲棒束區3個分區類型，在初始計算前進行如下參數設置：進口區黏性阻力系數初始指定為0，慣性阻力系數初始指定為0；出口區黏性阻力系數指定為0，慣性阻力系數指定為0，并在計算過程中保持不變；繞絲棒束區黏性阻力系數和慣性阻力系數通過經驗關系式獲得，并在計算過程中保持不變。

在穩態調試中，每隔n個迭代步獲取堆芯每盒組件的出口流量，并根據式(1)計算堆芯內各組件出口流量W的設計值與計算值之間的相對誤差εi。

(1)

式中，下標i為堆芯內組件序號，cal為計算值，des為設計值。

根據式(1)計算堆芯內各組件出口流量的設計值與計算值之間的相對誤差εi，通過式(2)調整堆芯各組件進口區多孔介質的黏性阻力系數和慣性阻力系數，直至各組件的流量與設計流量的誤差在允許范圍，認為穩態調試完成，可開展瞬態計算。

(2)

式中：1/α為黏性阻力系數；C2為慣性阻力系數；下標n為當前迭代步，m為上一迭代步，u為用戶指定修正的參數。

2) 物理模型

本文湍流模型采用Realizablek-ε模型、標準壁面函數、SIMPLE算法，梯度離散格式為Green-Gauss Node Based格式，壓力離散格式為二階格式，動量、湍動能、湍流耗散率和能量方程采用二階迎風格式，時間離散格式采用一階顯式格式，由于密度變化不大，浮升力作用采用Boussinesq假設。

3) 邊界條件和初始條件

計算時，進口采用流量進口邊界和進口溫度邊界，出口采用壓力邊界。其中進口溫度邊界采用圖7所示的瞬態溫度變化，該瞬態曲線根據實際反應堆實驗數據獲得。進口流量變化(歸一化值)如圖8所示，該曲線由實驗所得，該值乘以穩態流量即為瞬態流量入口邊界。堆芯各類組件的穩態功率分布如圖9所示，而瞬態下堆芯功率變化曲線(歸一化值)如圖8所示，該值乘以穩態功率即為瞬態功率變化值。需要指出的是，在瞬態計算過程假設堆芯各組件功率份額與穩態下保持一致。

圖7 冷管道進口瞬態溫度變化Fig.7 Transient temperature change at inlet of cold leg

圖8 歸一化堆芯功率和入口流量Fig.8 Normalized of core power and cold leg inlet mass flow rate

圖9 穩態堆芯功率分布Fig.9 Core power distribution under steady condition

2 計算結果與分析

2.1 穩態計算結果

每盒組件經過調整多孔介質參數后流量的計算值與設計值的相對誤差示于圖10，每盒組件的流量的相對誤差均小于±2%，且絕大多數組件流量的相對誤差小于±1%，因此可認為經過調整后的多孔介質較為合理，能滿足堆芯的流量分配需求，堆芯穩態計算結果是合理的。穩態下堆芯各類組件的出口溫度統計示于圖11，可看出，堆芯組件的最高出口溫度約為699 K，最低出口溫度約為602 K。

2.2 失流事故計算結果

1) 熱池流動換熱特性分析

反應堆堆芯出口溫度即熱管段溫度的模擬值與實驗值的對比示于圖12，可看出，鈉池巨大的鈉裝量使得在整個瞬態變化過程中溫度實驗值變化很小，前200 s內基本保持不變，而后700 s緩慢下降，且在整個瞬態過程僅有4 K的溫度下降。CFD模擬值較好地捕捉了反應堆出口溫度的變化趨勢。

圖10 組件流量計算值與設計值的相對誤差Fig.10 Relative error between calculated value and design value of assembly mass flow rate

圖11 穩態下堆芯組件出口溫度Fig.11 Outlet temperature of assembly under steady condition

圖12 熱管段溫度對比Fig.12 Comparison for hot leg temperatures

堆本體xz截面的溫度分布云圖示于圖13。從圖13可看出，初始階段，由于堆芯組件流速較高，致使在熱池中心區域形成1個高溫區域，而熱池外圍溫度相對較低，當堆芯流速降低后，中心高溫區域逐漸萎縮，從熱池頂部至熱池底部逐漸形成明顯的熱分層現象，特別是在熱池底部，由于流速低，從堆芯出來的冷流體可能積聚于此，直接導致此處的熱分層最為明顯。而對于冷池與其他腔室，整體溫度變化不大。

熱池xz截面不同時刻的速度矢量圖和速度云圖示于圖14。在初始階段，堆芯出口的流體向上流到液面處，在熱池中心區域形成1個高速區，而在中心高速區外，形成了2個大的漩渦，而后從熱管流出，同時伴隨著流速的下降，中心高速區域逐漸減小，特別是在后期，流速較低時，熱池上部的流體幾乎保持靜止，流體在幾乎與出口相同高度位置向外流出反應堆。

圖13 xz截面溫度云圖Fig.13 Temperature contour of xz cross section

圖14 熱池xz截面的速度矢量圖和速度云圖Fig.14 Vector and contour of velocity for xz cross section of hot pool

2) 堆芯流動換熱特性分析

具有溫度測點的組件的出口溫度模擬值與實驗值的對比示于圖15，該組件在堆芯中的位置如圖3所示。對于Row 6組件的溫度預測，實驗值除出現2個峰值溫度外，在100～150 s左右，溫度出現1個較為平緩的趨勢，而后在150 s后又開始下降，模擬值沒有體現該溫度變化趨勢，但考慮到Row 6組件在堆芯中的位置較特殊，周圍為氣體膨脹組件和非燃料組件，是堆內熱工和物理條件最為復雜的組件，從熱工角度來說，外側非燃料組件溫度較低，而內側為燃料組件，溫度較高，此處溫度梯度較大，盒間的導熱或對流換熱作用較強烈；從物理角度來說，瞬態下，由于氣體膨脹組件失去強迫循環，氣體膨脹組件液位下降，堆芯中子泄漏增加，可能導致該組件的功率份額分布與穩態分布差距較大，屬于物理熱工強耦合作用，單獨考慮熱工影響，可能不能準確模擬該過程。

圖15 Row 6 PIOTA組件出口溫度Fig.15 Outlet temperature of Row 6 PIOTA

堆芯活性出口截面的溫度分布示于圖16。可看出，堆芯最熱組件的位置在堆芯中不斷發生變化，剛開始位于堆芯中心區域，而后靠近燃料組件的外圍區域，最終呈現出中心低、中間環型區域高、外側低的馬鞍形溫度分布，而從堆芯y=0截面的溫度分布曲線可看出，剛開始時堆芯功率較高，y=0截面溫度分布差距較大，而伴隨瞬態進程，由于盒間的存在，溫度分布逐漸被展平。

不同時刻下盒間出口的流動特性示于圖17。在初始階段，從外側非燃料組件區域向下流入盒間，而從中心較高溫區的燃料區域流出盒間，而在后期，由于堆芯溫度分布產生了較大變化，使得盒間從外部和中心區域流入堆芯，而從外側燃料區域流出堆芯。通過對比盒間的矢量圖和堆芯xy截面的溫度分布云圖，可看出矢量圖流動方向基本與溫度保持一致，都是從低溫區流入盒間，從高溫區流出盒間，這是由于低溫時密度較大，而高溫時密度較小，在浮升力的作用下，呈現這種流動趨勢。

圖17 盒間出口截面流動特性Fig.17 Flow characteristic of outlet for inter wrapper

3 結論

本文采用CFD方法針對FFTF的LOSWOS#13無保護失流事故進行了900 s的瞬態模擬，并將模擬值與實驗值進行了對比，得出如下結論：

1) CFD方法能較好地捕捉鈉冷快堆全堆的熱池內熱分層等復雜的三維流動換熱特性，特別是對熱池巨大的熱容及其在瞬態下的熱惰性響應；

2) 堆芯最熱組件在瞬態過程中的位置不斷發生變化，堆芯盒間流的流動方向受盒內溫度影響較大，且盒間在堆芯流速低、功率小的情形下，能有效展平堆芯溫度。