基于CFD程序的物理熱工耦合計(jì)算不確定性分析

于 濤，雷洲陽，趙鵬程,*，錢冠華，李 捷

(1.南華大學(xué) 核科學(xué)技術(shù)學(xué)院，湖南 衡陽 421001；2.南華大學(xué) 湖南省數(shù)字化反應(yīng)堆工程技術(shù)研究中心，湖南 衡陽 421001)

物理熱工耦合方法中的不確定性量化是保障物理熱工耦合計(jì)算結(jié)果可靠性的重要基礎(chǔ)，是物理熱工耦合的重要研究方向之一。不確定性分析主要采用輸入不確定性傳播方法，研究輸入?yún)?shù)對(duì)模型帶來的不確定性[1]。將不確定性分析方法應(yīng)用到物理熱工耦合計(jì)算，能更精確地預(yù)估模型的參數(shù)范圍，為核反應(yīng)堆設(shè)計(jì)、安全分析和評(píng)審提供可靠的計(jì)算置信區(qū)間，具有重要的學(xué)術(shù)意義和實(shí)用價(jià)值[2-3]。

國(guó)際上，由美國(guó)橡樹嶺國(guó)家實(shí)驗(yàn)室(ORNL)牽頭的輕水堆先進(jìn)仿真聯(lián)盟(CASL)于2011年起采用DAKOTA軟件開展對(duì)多物理耦合的不確定性分析，取得了較為豐富的成果[4-9]；美國(guó)能源部(DOE)核能辦公室成立的核能先進(jìn)建模與仿真項(xiàng)目(NEAMS)于2017—2019年采用PROTEUS-Nek5000耦合工具箱SHARP開展了鈉冷快堆和鉛冷快堆安全特性的不確定性分析研究[10]；歐盟在“歐洲核能技術(shù)平臺(tái)可持續(xù)發(fā)展(SNETP)”戰(zhàn)略規(guī)劃中開展了NURESIM系列項(xiàng)目，在該項(xiàng)目中基于SALOME平臺(tái)進(jìn)行了多物理、多尺度耦合相關(guān)研究，并采用URANIE工具箱開展了相關(guān)不確定性分析工作[11-13]。近年來，國(guó)內(nèi)也進(jìn)行了一系列不確定性研究：沈昊等[14]采用敏感性系數(shù)法對(duì)棒束單相流動(dòng)傳熱計(jì)算流體力學(xué)(CFD)進(jìn)行了不確定性量化分析；潘昕懌等[15]采用混合法和隨機(jī)抽樣法研究了多群核數(shù)據(jù)不確定性對(duì)堆芯物理計(jì)算的影響；高新力等[16]基于SNAP-DAKOTA-RELAP5程序開展了壓水堆大破口事故的不確定性分析；張梟羽等[17]采用多項(xiàng)式混沌法進(jìn)行了壓水堆堆芯冷卻劑通道計(jì)算的不確定性分析；趙陽[18]基于DAKOTA/SALib程序開展了釷基熔鹽堆系統(tǒng)事故不確定性分析；曹志偉等[19]采用確定論統(tǒng)計(jì)法(GSM)對(duì)CPR1000核電廠大破口事故進(jìn)行了不確定性量化分析；郭炯等[20]基于VSOP-UAM程序開展了高溫氣冷堆的不確定性分析；劉勇等[21]采用兩步法的敏感性計(jì)算策略，針對(duì)快堆BN-600基準(zhǔn)例題的有效增殖因數(shù)進(jìn)行了不確定性量化分析；鄧程程等[22]基于RELAP5最佳估算程序?qū)ο冗M(jìn)熱工水力試驗(yàn)(ACME)臺(tái)架開展了小破口失水事故的不確定性和敏感性分析。一方面，國(guó)內(nèi)在堆物理、系統(tǒng)程序分析不確定性量化方面進(jìn)行了大量工作，但有關(guān)物理熱工耦合不確定性分析工作開展較少；另一方面，第4代反應(yīng)堆中鈉冷快堆以及鉛基快堆存在極為復(fù)雜的物理熱工耦合現(xiàn)象，應(yīng)用不確定性分析可更精確地評(píng)估反應(yīng)堆安全特性，因此有必要進(jìn)行深入研究。

本文基于CFD程序FLUENT的用戶自定義函數(shù)(UDF)，耦合中子動(dòng)力學(xué)計(jì)算模型、燃料棒熱傳導(dǎo)計(jì)算模型、不確定性分析程序SIMLAB，開發(fā)物理熱工耦合計(jì)算不確定性分析平臺(tái)CFD/PFS。進(jìn)而基于CFD/PFS開展小型自然循環(huán)鉛冷快堆SNCLFR-10的無保護(hù)超功率(UTOP)事故不確定性分析，研究燃料多普勒系數(shù)等不確定性輸入?yún)?shù)對(duì)功率峰值等目標(biāo)參數(shù)帶來的不確定性，并評(píng)估反應(yīng)性反饋系數(shù)對(duì)目標(biāo)參數(shù)的影響程度。

1 不確定性分析CFD/PFS平臺(tái)開發(fā)

1.1 點(diǎn)堆中子動(dòng)力學(xué)模型

點(diǎn)堆中子動(dòng)力學(xué)模型(PKM)主要用于計(jì)算反應(yīng)堆內(nèi)中子密度的動(dòng)態(tài)變化，進(jìn)而更新堆芯功率，為冷卻劑和慢化劑熱工水力計(jì)算提供熱源[23]。本工作采用6組緩發(fā)中子的點(diǎn)堆中子動(dòng)力學(xué)模型：

(1)

(2)

其中：N(t)為中子密度；ρ(t)為反應(yīng)性；Λ為中子代時(shí)間；U為緩發(fā)中子總份額；t為時(shí)間；Ui為第i組緩發(fā)中子份額；λi為第i組緩發(fā)中子先驅(qū)核衰變常量；Ci(t)為第i組緩發(fā)中子先驅(qū)核濃度；q為外加中子源源強(qiáng)。采用二階泰勒多項(xiàng)式積分法求解點(diǎn)堆中子動(dòng)力學(xué)方程，將中子密度N(t)使用全隱式二階泰勒近似展開，以迭代的方式求解不同時(shí)刻的中子密度。

點(diǎn)堆中子動(dòng)力學(xué)模型中，中子密度主要由總反應(yīng)性控制，總反應(yīng)性包括外部引入反應(yīng)性和內(nèi)部反應(yīng)性反饋[24]。反饋反應(yīng)性主要包括由燃料溫度引起的燃料多普勒反應(yīng)性、材料密度變化引起的軸向膨脹反應(yīng)性和徑向膨脹反應(yīng)性以及冷卻劑密度反應(yīng)性?？偡磻?yīng)性為：

αc,r(Tc(t)-Tc(0))+αc(Tc(t)-Tc(0))

(3)

其中：KD為多普勒系數(shù)；αc,a為堆芯軸向膨脹反饋系數(shù)；αc,r為堆芯徑向膨脹反饋系數(shù)；αc為冷卻劑溫度反饋系數(shù)；Tf(t)和Tc(t)分別為t時(shí)刻的燃料、冷卻劑平均溫度。

1.2 燃料棒熱傳導(dǎo)模型

燃料棒熱傳導(dǎo)模型(FTM)包括穩(wěn)態(tài)和瞬態(tài)熱傳導(dǎo)模型，用于求解燃料棒芯塊、氣隙、包殼的溫度，將包殼表面熱流密度傳遞給CFD程序，并反饋PKM燃料芯塊和包殼的溫度。

1) 穩(wěn)態(tài)熱傳導(dǎo)模型

芯塊熱傳導(dǎo)：

(4)

氣隙熱傳導(dǎo)：

(5)

包殼內(nèi)部熱傳導(dǎo)：

(6)

包殼與冷卻劑間的對(duì)流換熱：

Tcu-Tw=ql/(2π·Rcu·h)

(7)

換熱系數(shù)h的計(jì)算關(guān)系式：

(8)

積分熱導(dǎo)率[25]：

(9)

其中：Ru為芯塊半徑；Rci為包殼內(nèi)半徑；Rcu為包殼外半徑；T0和Tu分別為芯塊中心、外表面溫度；Tci和Tcu分別為包殼內(nèi)、外表面溫度；Tw為冷卻劑溫度；ku、kg、kz、kc依次為芯塊、間隙、包殼、冷卻劑的熱導(dǎo)率；h為包殼與冷卻劑間的換熱系數(shù)；ql為線功率密度；Qv為體積釋熱率；Nu為努塞爾數(shù)；De為當(dāng)量直徑；k(T)為對(duì)應(yīng)材料的熱導(dǎo)率公式；k為材料的平均熱導(dǎo)率；T為材料溫度；Ta、Tb為相鄰兩節(jié)點(diǎn)溫度。采用積分熱導(dǎo)率法求解燃料棒穩(wěn)態(tài)溫度，根據(jù)邊界條件由外向內(nèi)迭代求解燃料包殼、芯塊的溫度分布。

2) 瞬態(tài)熱傳導(dǎo)模型

根據(jù)能量守恒，燃料棒瞬態(tài)熱傳導(dǎo)模型的控制方程[26]為：

(10)

忽略燃料棒軸向傳熱，則：

(11)

(12)

其中：下標(biāo)m代表材料，包括芯塊、氣隙、包殼材料；ρm為材料密度；cp,m為材料比定壓熱容；r為半徑；k為材料熱導(dǎo)率；i為徑向所在的節(jié)點(diǎn)數(shù)。采用有限差分方法求解燃料棒瞬態(tài)導(dǎo)熱模型，沿徑向?qū)⑿緣K和包殼劃分多個(gè)節(jié)點(diǎn)，迭代求解燃料棒瞬態(tài)溫度分布。

1.3 計(jì)算流體力學(xué)模型

CFD模型用于計(jì)算冷卻劑溫度、速度等參數(shù)，將冷卻劑溫度、速度傳遞給FTM，并反饋PKM冷卻劑的溫度。CFD模型中的控制方程包括質(zhì)量守恒方程、動(dòng)量守恒方程和能量守恒方程，表達(dá)式如下。

質(zhì)量守恒方程：

(13)

動(dòng)量守恒方程：

(14)

能量守恒方程：

(15)

實(shí)際反應(yīng)堆的幾何結(jié)構(gòu)非常復(fù)雜，為降低建模難度、節(jié)省計(jì)算時(shí)間，采用多孔介質(zhì)模型分析堆芯和熱交換器模型，將冷卻劑流經(jīng)這些區(qū)域的壓降等效為動(dòng)量守恒方程中的附加動(dòng)量源項(xiàng)[27]，表達(dá)式為：

(16)

其中：u為速度矢量；t為時(shí)間；ρ為液體密度；ui為速度場(chǎng)分量；μ為湍動(dòng)能黏度；p為液體壓力；xi為方向坐標(biāo)；Si為外加動(dòng)力源項(xiàng)；k為流體的熱導(dǎo)率；cp為比定壓熱容；ST為體積熱源項(xiàng)；vi為局部流速；μ(T)為隨溫度變化的流體動(dòng)力黏度；ρ(T)為隨溫度變化的流體密度；Di為黏性壓降系數(shù)；Ci為慣性壓降系數(shù)。采用有限體積法(FVM)建立離散控制方程組，求解離散方程時(shí)選取SIMPLE算法計(jì)算流體的溫度與速度[28]。

1.4 不確定性分析程序SIMLAB

SIMLAB是由歐洲聯(lián)合研究中心(JRC)研發(fā)的一款不確定性開源程序，具有可視化操作界面，包括統(tǒng)計(jì)預(yù)處理模塊、模型執(zhí)行模塊、統(tǒng)計(jì)后處理模塊，主要用于不確定性分析和全局敏感性分析[29]。SIMLAB程序在管理經(jīng)濟(jì)可行性評(píng)估以及模糊推理系統(tǒng)等領(lǐng)域已經(jīng)開展了不確定性和敏感性分析[30-31]。

SIMLAB程序的統(tǒng)計(jì)預(yù)處理模塊提供參數(shù)分布類型及抽樣方法，用戶根據(jù)模型要求選擇不同的參數(shù)分布及抽樣方法。SIMLAB程序中參數(shù)分布類型包括正態(tài)分布、均勻分布、韋伯分布、指數(shù)分布、γ分布、β分布等，也可以自行創(chuàng)建函數(shù)關(guān)系，且能可視化參數(shù)分布。抽樣方法主要包括FAST法、Sobol法、拉丁超立方方法、隨機(jī)抽樣法等。

使用非參數(shù)統(tǒng)計(jì)法(Wilks公式)確定抽樣次數(shù)，Wilks公式根據(jù)分布區(qū)間的置信度b及輸出參數(shù)的概率a確定所需最小抽樣數(shù)目n[32]。Wilks公式為：

1-an≥b

(17)

(1-an)-n(1-a)an-1≥b

(18)

其中式(17)適用于單側(cè)限值分布，式(18)適用于雙側(cè)限值分布。

SIMLAB程序的統(tǒng)計(jì)后處理模塊用于不確定性分析及敏感性評(píng)估。不確定性分析能可視化參數(shù)的概率分布、累積分布、逆累積分布，給出參數(shù)的均值、方差、標(biāo)準(zhǔn)差、上下限值等統(tǒng)計(jì)結(jié)果；敏感性評(píng)估給出輸入?yún)?shù)與輸出參數(shù)的相關(guān)系數(shù)或敏感性指數(shù)等。

1.5 物理熱工耦合不確定性分析平臺(tái)建立

熱工計(jì)算包括計(jì)算流體力學(xué)模型和燃料棒熱傳導(dǎo)模型，物理計(jì)算為點(diǎn)堆中子動(dòng)力學(xué)模型，3個(gè)模型的耦合方法如圖1所示。PKM進(jìn)行堆芯中子學(xué)計(jì)算，對(duì)FTM提供堆芯功率；FTM進(jìn)行瞬態(tài)熱傳導(dǎo)計(jì)算，對(duì)CFD提供熱流密度，對(duì)PKM提供燃料溫度；CFD進(jìn)行流動(dòng)換熱計(jì)算，對(duì)FTM提供冷卻劑的溫度和速度。

圖1 物理熱工耦合方法

物理熱工耦合程序計(jì)算流程如圖2所示，首先讀取初始、邊界條件，調(diào)用CFD模塊計(jì)算穩(wěn)態(tài)的冷卻劑溫度、速度，將CFD模塊計(jì)算結(jié)果傳遞給FTM計(jì)算燃料棒溫度分布，當(dāng)前后兩次計(jì)算的各節(jié)點(diǎn)溫度分布小于0.1 K時(shí)，輸出穩(wěn)態(tài)計(jì)算值；讀取穩(wěn)態(tài)計(jì)算結(jié)果作為瞬態(tài)初始參數(shù)，將初始化參數(shù)傳遞給PKM計(jì)算堆芯功率；FTM通過堆芯功率計(jì)算燃料棒溫度分布和表面熱流密度，CFD模塊根據(jù)表面熱流密度計(jì)算冷卻劑溫度、速度分布；將燃料棒、冷卻劑溫度分布反饋給PKM用于計(jì)算反應(yīng)性反饋，并將冷卻劑溫度、速度反饋給FTM。當(dāng)計(jì)算時(shí)間達(dá)到耦合時(shí)間步長(zhǎng)，物理與熱工模塊進(jìn)行數(shù)據(jù)傳遞，之后進(jìn)行下一時(shí)間步的計(jì)算，直至FTM計(jì)算出的燃料、包殼溫度幾乎不變，則認(rèn)為瞬態(tài)事故區(qū)域穩(wěn)定，停止計(jì)算。

圖2 物理熱工耦合程序計(jì)算流程

通過SIMLAB程序中模型執(zhí)行模塊的外部接口，實(shí)現(xiàn)物理熱工耦合程序與SIMLAB程序的耦合計(jì)算，開發(fā)了物理熱工耦合不確定性分析平臺(tái)CFD/PFS。耦合方法如圖3所示：通過SIMLAB程序?qū)︸詈铣绦蛑械牟淮_定性輸入?yún)?shù)進(jìn)行隨機(jī)抽樣，并生成樣本文件；將樣本文件傳遞給物理熱工耦合程序，經(jīng)計(jì)算獲得輸出文件；通過執(zhí)行模塊外部接口，將輸出文件傳遞回SIMLAB程序，并應(yīng)用SIMLAB程序統(tǒng)計(jì)分析樣本參數(shù)與輸出參數(shù)的相關(guān)系數(shù)，對(duì)參數(shù)進(jìn)行不確定性量化分析。

基于CFD/PFS平臺(tái)的不確定性計(jì)算流程如圖4所示，計(jì)算流程如下：1) 識(shí)別重要的不確定性輸入?yún)?shù)；2) 在SIMLAB程序中定義重要輸入?yún)?shù)的取值范圍及概率分布，選取抽樣方法，獲得不同的樣本組合；3) 將每一組樣本組合傳遞給物理熱工耦合程序進(jìn)行計(jì)算獲得響應(yīng)值；4) 將響應(yīng)值傳遞給SIMLAB程序，利用SIMLAB程序?qū)Y(jié)果進(jìn)行不確定性量化，計(jì)算輸入?yún)?shù)與響應(yīng)值的相關(guān)系數(shù)。

圖3 SIMLAB與物理熱工耦合程序的耦合方法

圖4 不確定性方法的計(jì)算流程

2 物理耦合程序驗(yàn)證

2.1 分析模型

小型自然循環(huán)鉛基快堆(SNCLFR-10)是由中國(guó)科學(xué)院設(shè)計(jì)的功率為10 MW的鉛鉍合金冷卻的池式反應(yīng)堆，燃料采用富集度為19.75%的UO2，燃料包殼采用316L不銹鋼，一回路冷卻劑采用液態(tài)鉛鉍合金共金體[33]。SNCLFR-10的結(jié)構(gòu)示意圖示于圖5，堆芯主要設(shè)計(jì)參數(shù)參見文獻(xiàn)[34]。采用國(guó)際著名快堆多物理耦合分析程序SIMMER-Ⅲ，基于SNCLFR-10無保護(hù)超功率(UTOP)事故對(duì)CFD/PFS平臺(tái)物理熱工耦合程序開展Code-to-Code驗(yàn)證。

圖5 SNCLFR-10的結(jié)構(gòu)

2.2 耦合程序驗(yàn)證

瞬態(tài)超功率事故是反應(yīng)堆設(shè)計(jì)必須考慮的典型事故瞬態(tài)工況，是指向堆內(nèi)突然引入一個(gè)外部反應(yīng)性，導(dǎo)致功率急劇上升的事故。本文將反應(yīng)堆初始運(yùn)行條件設(shè)為穩(wěn)態(tài)，在2 s內(nèi)線性引入1 $反應(yīng)性，整個(gè)瞬態(tài)過程無停堆動(dòng)作。

UTOP事故中反應(yīng)性的計(jì)算值和功率變化如圖6所示。反應(yīng)性在2 s內(nèi)激增，但由于燃料和冷卻劑的負(fù)反饋?zhàn)饔?，反?yīng)性只能提高到0.88 $左右，之后由于燃料和冷卻劑的負(fù)反饋增大，反應(yīng)性逐漸降低，反應(yīng)性最后穩(wěn)定在臨界狀態(tài)。反應(yīng)堆功率在2 s內(nèi)升高30倍左右，由于負(fù)反饋?zhàn)饔茫β手饾u降低，最后穩(wěn)定在58 MW左右。SIMMER-Ⅲ的功率最后穩(wěn)定在60 MW左右，耦合程序與SIMMER-Ⅲ的相對(duì)誤差為3.33%。

圖7示出燃料最高溫度和包殼最高溫度隨時(shí)間的變化。隨著功率的增大，燃料最高溫度和包殼最高溫度逐漸升高，之后隨著負(fù)反應(yīng)性的增大，堆芯功率降低，燃料和包殼的最高溫度會(huì)隨之降低。相對(duì)于功率的變化，燃料和包殼的最高溫度變化幅度較低，其變化在時(shí)間上也比功率延遲。耦合程序和SIMMER-Ⅲ計(jì)算的燃料最高溫度最后分別穩(wěn)定在2 500 K和2 450 K左右，兩者的相對(duì)誤差為2.04%；包殼最高溫度分別穩(wěn)定在1 098 K和1 080 K左右，兩者的相對(duì)誤差為1.7%。

圖6 反應(yīng)性和堆芯功率隨時(shí)間的變化

圖7 燃料最高溫度和包殼最高溫度隨時(shí)間的變化

對(duì)比物理熱工耦合程序與SIMMER-Ⅲ的計(jì)算值，兩者的誤差較小、符合度較高，說明CFD/PFS平臺(tái)計(jì)算結(jié)果具有良好的準(zhǔn)確性和可靠性，可用于開展反應(yīng)堆物理熱工耦合不確定性量化分析研究。

3 UTOP事故不確定性分析

3.1 不確定性參數(shù)的選取

反應(yīng)堆中，熱工水力反饋參數(shù)會(huì)通過影響相關(guān)核素的核密度或微觀截面，從而直接影響組件等效均勻化宏觀截面參數(shù)，因此將燃料多普勒系數(shù)、堆芯軸向膨脹系數(shù)、包殼徑向膨脹系數(shù)和冷卻劑溫度系數(shù)作為不確定性輸入?yún)?shù)[35]。為減少抽樣取值的偶然性，假設(shè)所有輸入?yún)?shù)為均勻分布；將輸入?yún)?shù)的相對(duì)不確定度假設(shè)為5%。輸入?yún)?shù)的不確定性分布列于表1。

表1 輸入?yún)?shù)的不確定性分布

在SNCFLR-10的UTOP事故分析中，將表征反應(yīng)堆運(yùn)行狀態(tài)及安全特性的參數(shù)，包括總反應(yīng)性峰值、功率峰值、燃料及包殼最高溫度，選取為不確定性輸出參數(shù)。

3.2 不確定性量化

選取拉丁超立方方法對(duì)輸入?yún)?shù)抽樣，根據(jù)式(18)，滿足兩個(gè)95%雙側(cè)限值分布的最小樣本量為93。抽取100組工況作為輸入樣本，為直觀檢測(cè)輸入?yún)?shù)的樣本量，將輸入?yún)?shù)進(jìn)行歸一化處理[36]。假設(shè)第i個(gè)參數(shù)的上限值為Uij，下限值為L(zhǎng)ij，第j個(gè)樣本的參數(shù)值為Kij，則歸一化結(jié)果Xij為：

(19)

Xij的取值在0～1之間，靠近0的樣本接近下限值，靠近1的樣本接近上限值。

圖8示出拉丁超立方抽樣樣本歸一化分布。由圖8可見，歸一化后的參數(shù)服從均勻分布，因此輸入?yún)?shù)也服從均勻分布，驗(yàn)證了輸入?yún)?shù)抽樣的合理性。

圖8 拉丁超立方抽樣樣本歸一化分布

圖9示出總反應(yīng)性峰值、功率峰值、包殼最高溫度、燃料最高溫度的分布。由圖9可知：總反應(yīng)性峰值的上限值為0.889 13 $，下限值為0.883 89 $，名義值為0.886 63 $，名義值與限值的最大相對(duì)偏差為3.09%；功率峰值的上限值為322.2 W，下限值為298.04 W，名義值為310.29 MW，名義值與限值的最大相對(duì)偏差為3.948%；燃料最高溫度的上限值為2 746.966 K，下限值為2 560.997 K，名義值為2 656.86 K，名義值與限值的最大相對(duì)偏差為3.608%；包殼最高溫度的上限值為1 221.189 K，下限值為1 169.641 K，名義值為1 196.409 K，名義值與限值的最大相對(duì)偏差為2.237%。各目標(biāo)參數(shù)的名義值處于合理的結(jié)果不確定性范圍之內(nèi)，表明選取的輸入?yún)?shù)分布是合理的；名義值與限值的相對(duì)偏差均不超過3.95%，表明計(jì)算的目標(biāo)參數(shù)精度較高。

3.3 敏感性分析

輸入?yún)?shù)之間的相關(guān)性可能影響輸入和輸出的相關(guān)系數(shù)，利用偏相關(guān)系數(shù)(PCC)對(duì)計(jì)算結(jié)果進(jìn)行分析。偏相關(guān)系數(shù)表示消除其他影響因素后輸入?yún)?shù)對(duì)響應(yīng)值的影響程度，取值在-1～1之間。PCC為0時(shí)，輸入變量與響應(yīng)值不相關(guān)；PCC為正呈正相關(guān)，為負(fù)呈負(fù)相關(guān)；PCC絕對(duì)值靠近1時(shí)，相關(guān)性最強(qiáng)。假設(shè)給定隨機(jī)變量X1和X2作為輸入，輸出變量Y，在消除X2所產(chǎn)生的影響后，PCC計(jì)算表達(dá)式為：

圖9 總反應(yīng)性峰值、功率峰值、燃料最高溫度、包殼最高溫度的分布

(20)

其中：RX1Y|X2為偏相關(guān)系數(shù)；RX1Y為X1與Y的簡(jiǎn)單相關(guān)系數(shù)；RX1X2為X1與X2的簡(jiǎn)單相關(guān)系數(shù)變量間的相關(guān)系數(shù)；RX2Y為X2與Y的簡(jiǎn)單相關(guān)系數(shù)。

圖10示出輸入?yún)?shù)與響應(yīng)值的偏相關(guān)系數(shù)。由圖10可見，輸入?yún)?shù)與響應(yīng)值均為正相關(guān)，其中，冷卻劑溫度系數(shù)與總反應(yīng)性峰值、功率峰值、燃料最高溫度及包殼最高溫度的偏相關(guān)系數(shù)均小于0.16，屬于極弱相關(guān)；包殼徑向膨脹系數(shù)、堆芯軸向膨脹系數(shù)與各響應(yīng)值的偏相關(guān)系數(shù)在0.4～0.8之間，相關(guān)性較強(qiáng)；燃料多普勒系數(shù)與各響應(yīng)值的偏相關(guān)系數(shù)大于0.94，屬于極強(qiáng)相關(guān)。

圖10 輸入?yún)?shù)與響應(yīng)值的偏相關(guān)系數(shù)

4 結(jié)論

針對(duì)物理熱工耦合計(jì)算進(jìn)行不確定性分析，能得到更精確的關(guān)鍵參數(shù)分布信息，可進(jìn)一步提高反應(yīng)堆的安全性。本文介紹了基于CFD程序FLUENT開發(fā)的不確定性分析平臺(tái)CFD/PFS，并采用國(guó)際著名快堆多物理耦合分析程序SIMMER-Ⅲ對(duì)CFD/PFS平臺(tái)開展對(duì)比驗(yàn)證，進(jìn)而基于CFD/PFS平臺(tái)對(duì)SNCLFR-10的UTOP事故開展了不確定性分析，并研究了不確定性參數(shù)的敏感性，所得結(jié)論如下。

1) CFD/PFS平臺(tái)和SIMMER-Ⅲ程序的事故分析結(jié)果吻合良好，所開發(fā)的CFD/PFS平臺(tái)計(jì)算結(jié)果具有較好的準(zhǔn)確性和可靠性，可用于物理熱工耦合不確定性量化分析。

2) 通過不確定性分析計(jì)算，得到了SNCLFR-10的UTOP事故中總反應(yīng)性峰值、功率峰值、燃料最高溫度及包殼最高溫度等瞬態(tài)安全參數(shù)的95/95雙側(cè)容忍限值，且瞬態(tài)安全參數(shù)的名義值均處于雙側(cè)容忍限值內(nèi)，表明對(duì)UTOP事故的不確定性分析具有合理性。量化分析中名義值與限值的相對(duì)偏差小于3.95%，表明獲得的瞬態(tài)安全參數(shù)不確定性較小，計(jì)算結(jié)果精度較高。

3) 對(duì)于SNCLFR-10的UTOP事故敏感性分析，通過對(duì)4組反應(yīng)性反饋系數(shù)的偏相關(guān)系數(shù)比較得知，堆芯軸向膨脹系數(shù)、包殼徑向膨脹系數(shù)及燃料多普勒系數(shù)對(duì)瞬態(tài)安全參數(shù)的影響較為顯著。其中燃料多普勒系數(shù)對(duì)結(jié)果的影響極為強(qiáng)烈，對(duì)反應(yīng)堆安全影響最大。