COBRA-Ⅳ對8×8棒束計算的不確定性分析

杜 蕓，劉曉晶，程 旭

(上海交通大學 核能科學與工程學院，上海 200240)

1988年，美國NRC提出了最佳估算模型的安全分析程序，針對最佳估算程序的不確定性分析方法也應運而生。不確定性分析成為如今進行核能熱工安全分析必不可少的工作。然而，由于核能領域的不確定性分析方法最初是針對系統程序展開的，針對子通道程序開發的不確定性分析方法還較少。子通道分析程序是計算反應堆堆芯熱工水力過程現象的通用工具，子通道程序的模擬計算同樣存在不確定性，對其進行不確定性分析十分必要。

目前針對系統程序開發的不確定性分析方法[1]，從跟蹤不確定性的方法上可分為輸入參數不確定性的傳播和輸出結果的誤差傳播兩類。輸入參數不確定性的傳播可通過如下方式得到：首先認定并給出不確定性參數的范圍和分布，然后通過改變這些輸入參數進行計算。輸出結果的誤差傳播可通過計算結果與實驗數據的比較直接得到。本工作采用輸入參數不確定性傳播法[2]對COBRA-Ⅳ程序的計算進行不確定性分析。

1 計算對象

圖1 棒束子通道分布示意圖

本文研究的對象為BFBT[3]眾多實驗中的一個。BFBT是由美國NRC與日本金融、貿易和工業部一起核準，最后被經濟合作組織(OECD)認可的一國際性工程。本實驗為模擬沸水堆燃料棒束，建立一在高壓、高溫條件下垂直的8×8棒束，棒束橫截面和各子通道分布示意圖示于圖1。組件盒內裝有60根燃料棒，呈8×8方式排列，組件中間有一直徑為34.0 mm的不加熱的水棒，其中無流體流動。棒束的軸向功率均勻分布，徑向功率非均勻分布，徑向功率分布示于圖2。定位格架和棒束的一些相關參數列于表1。

圖2 徑向相對功率分布

表1 組件的幾何參數

本工作模擬計算的是該實驗中質量含氣率最高的1組穩態工況，該工況的參數列于表2。

表2 子通道計算工況

2 空泡份額計算值的影響因素

利用子通道分析程序COBRA-Ⅳ對以上實驗對象進行分析計算。如圖1所示，將冷卻劑流通橫截面劃分為80個子通道。為研究影響空泡份額計算值的因素，采用較簡單的模型對實驗進行模擬計算，作為基準算例進行比較。本工作從邊界條件和計算模型兩方面分析影響空泡份額計算值的因素。

2.1 邊界條件對空泡份額的影響

雖然實驗的邊界條件由實驗設備控制且由儀器測量，理論上符合工況的設定，但測量儀器存在誤差[3]。采用單一變量原則，即每次計算時只有1個邊界條件是變量，其余仍為工況設計值，以此逐一分析邊界條件的偏差對結果的影響。

為分析每個邊界條件的偏差對結果的影響并方便比較其對結果的影響，假設每個條件均存在同樣的相對偏差，綜合參考實驗設備對邊界條件的測量誤差[3]，將該相對偏差定為±1%。相對偏差的定義為：

×100%

(1)

其中：d為相對偏差；fa為實際值；fb為基準值。

分別計算邊界條件發生偏差時對子通道空泡份額計算結果造成的影響，得到相對變化量較大的子通道為31、33、48、50。這些子通道均是水棒周圍的子通道。

評估邊界條件的偏差對計算結果造成偏差的平均效應：

εre-d(i)×100%

(2)

其中：εre-d為相對變化量；i為子通道的通道編號。

分別評估邊界條件的偏差對空泡份額計算結果造成偏差的平均效應，結果列于表3。

表3 邊界條件對空泡份額的影響

由表3計算結果可見，對空泡份額影響最大的邊界條件是入口流體的焓，其次，出口壓力對空泡份額也有一定影響，而流量及熱流密度對空泡份額的影響較小。雖然邊界條件對空泡份額的影響不大，但邊界條件的不確定性是必然存在的，所以需考慮其對結果不確定性的影響。

2.2 計算模型對空泡份額的影響

程序中的模型及模型參數的選擇也會對計算結果產生影響。對COBRA-Ⅳ中涉及熱工水力計算的主要物理模型及模型參數進行研究。由于模型眾多，本文僅以空泡份額模型為例。

COBRA-Ⅳ中可供選擇的空泡份額模型[4]有以下幾種。

均相模型：

α=χvG/[χvG+(1-χ)vL]

(3)

Modified Armand模型：

α=(0.833+0.167χ)χvG/((1-χ)vL+χvG)

(4)

Chexal-Lellouche模型：

α=jG/[C0(jG+jL)+vgj]

(5)

滑移模型(滑速比可設置)：

α=χvG/[(1-χ)vLS+χvG]

(6)

其中：α為空泡份額；vG為氣體比體積；vL為液體比體積；χ為質量含氣率；jG為氣相表觀速度；jL為液相表觀速度；vgj為漂移速度；S為滑速比；C0為表示兩種速度關系的系數[4]。

基準計算采用均相模型，其余3種模型的計算結果與基準計算進行對比，其中，滑移模型中的滑速比是可設置的?；谧钚幽芰骷僭O，理論滑速比為2.7，因此，本文選擇的滑速比為2.0～3.0。

圖3 不同模型的空泡份額計算結果比較

不同模型空泡份額的計算結果比較示于圖3。圖3中，H、MA和CL分別代表均相模型、Modified Armand模型和Chexal-Lellouche模型，S2.0、S2.5和S3.0分別代表滑速比為2.0、2.5和3.0的滑移模型。由圖3可見：用不同模型計算的子通道出口空泡份額與基準計算結果的趨勢一致，只是大小有異；均相模型的計算結果最大，這一點由式(3)不難得出。同時，滑移模型的計算結果受滑速比的影響非常大，由于滑速比定義了流體中氣相與液相流動速度的關系，由式(6)可知該值在計算空泡份額時起重要作用。

圖3結果表明，空泡份額模型的選擇對計算結果的影響相當大。在本文計算范圍內，平均相對變化量高達20%，這比邊界條件對計算結果的影響大得多。可見，選擇適合工況的空泡份額模型十分重要。

對模型及模型參數依次進行研究，綜合所有模型對子通道空泡份額計算值造成的影響，將不同模型的計算結果與基準計算結果進行比較。取80個子通道中相對變化量的最大值的絕對值，結果列于表4。

由表4可見，與空泡份額的計算直接相關的空泡份額模型對計算結果的影響最大，空泡漂移流模型及修正種類的選擇對計算結果的影響也較大。交混系數對計算結果的影響較小，但該影響無法避免。

表4 各模型對子通道出口空泡份額計算值的影響比較

綜上所述，在確定采用滑移模型作為空泡份額模型后，確定6個參數為不確定性輸入參數。其中，4個參數為邊界條件，包括出口壓力、入口焓、入口質量流密度和熱流密度，2個參數為模型參數，包括滑速比和交混系數。邊界條件的不確定性范圍由儀器的測量精度決定，但不確定性的分布無從知道，所以均按照均勻分布來處理。本文考慮的輸入參數的不確定性范圍及分布列于表5。

表5 輸入參數的不確定性范圍及分布

3 不確定性分析

3.1 數學模型及原理

對于給定的計算容忍限，需要確定取樣數目[5]。文獻[6]提出對1個量進行雙邊容忍限(式(7))及單邊容忍限(式(8))的計算公式：

β=1-αN-N(1-α)αN-1

(7)

β=1-αN

(8)

其中：β為置信度；α為總體空間在兩個限值間的份額；N為最小采樣次數。

式(7)、(8)即為Wilks公式。

3.2 不確定性分析程序

確定了要求的容忍限，根據Wilks公式可得到需要的計算次數。該計算次數與不確定性輸入參數的個數、不確定性輸入參數的范圍及分布均無關，所以該方法適用性很廣。在滿足計算次數要求的同時，采樣方式須根據數學原理符合計算范圍以及分布，且隨機產生。這樣計算的結果經適當處理才能滿足要求的容忍限。

本工作編寫了子通道程序的不確定性分析程序，該程序的計算流程示于圖4。采用排序的方式對N組計算結果進行處理，將每個子通道的計算結果按升序排序，其中，排在中間的值作為該子通道的預測結果，排在第1位的值為滿足容忍限的不確定性下限，排在最后1位的值為不確定性的上限。這樣，便完成了對該子通道程序計算的不確定性分析。

圖4 不確定性分析程序的流程

3.3 程序驗證

運用本工作編寫的不確定性分析程序產生93組程序輸入參數。由于SUSA方法也使用了Wilks公式的原理，將本文不確定性分析程序產生的數據與SUSA方法的數據進行比較，驗證本方法的可行性。

1) 驗證各不確定性輸入參數是否符合分布規律

以入口焓采樣為例，用已排序的數據進行對比，結果示于圖5。由圖5可見，入口焓的采樣值分布近似呈線性，且采樣點疏密較一致，重合率較高，說明兩者采樣均是均勻分布。因此，本文的采樣能基本符合參數分布的特點。同樣，對于其余5個不確定性輸入參數，本文采樣與SUSA方法的也吻合較好，且符合每個參數分布的特點。

圖5 入口焓的采樣排序對比

2) 各參數組合比較

為比較6個輸入參數的組合方式，將每一種參數的93個采樣值劃分為大于基準值和小于等于基準值兩種，小于等于基準值的為0類，大于基準值的為1類。對程序產生的每一組參數組合按入口焓、入口流量、出口壓力、熱流密度、滑速比和交混系數的順序編號，表示每組參數組合的特征。如第1組參數的編號為011101，說明第1組參數的特征為：入口焓小于等于基準值，入口流量大于基準值，出口壓力大于基準值，熱流密度大于基準值，滑速比小于等于基準值，交混系數大于基準值。將編號作為一二進制數，便可轉換為一十進制數，共有64種組合方式，分別用0到63表示。這樣，第1組參數的特征可用29來表示。將SUSA和本文方法對于不確定性參數的組合用數字表示，并進行比較，結果示于圖6。

由圖6可知：SUSA方法覆蓋的組合類型為51種，本文方法覆蓋的組合類型為50種，覆蓋率均約為80%，二者均未覆蓋所有組合；本文的組合類型與SUSA方法的組合類型相比，相似率為80%。

綜上可見，本文采樣方法與組合方式是合格的，程序計算出的結果經順序處理，能滿足容忍限的要求。

圖6 不同種類組合數目的比較

3.4 不確定性分析結果

由于該實驗屬高空泡份額的實驗，滑移模型較適合計算這類空泡份額，因此選定滑移模型作為本計算的空泡份額模型。其余模型則根據棒束實驗，選擇較適合實驗的模型即可。將雙邊容忍限定為(95%，95%)，采用自行編寫的不確定性分析程序，得到最終的計算結果，即每個子通道的出口空泡份額預測值y，滿足(95%，95%)的不確定性上限y95/95上和不確定性下限y95/95下。

圖7示出子通道空泡份額的計算預測值與實驗值的對比。圖8示出空泡份額預測值與實驗值的相對誤差。由圖7、8可知，熱流密度較高的棒束周圍的子通道的空泡份額較高，說明子通道程序對空泡份額的預測符合實驗的規律以及趨勢。預測值普遍比實驗值小，平均相對誤差為10%，最大相對誤差在20%以內，均在可接受范圍內。與實驗值相對誤差較大的是受熱不均勻、不對稱的子通道。

圖7 子通道空泡份額預測值與實驗值的對比

對于堆芯棒束的熱工水力問題，最關心的是通道內冷卻劑的換熱能力能否得以保證，冷卻劑能否及時有效地將棒束產生的熱量帶走。空泡的存在會削弱冷卻劑的換熱能力，且空泡份額過大時有可能發生傳熱惡化現象，導致包殼溫度過高，堆芯安全受到威脅。觀察棒束的實驗結果以及計算結果(圖8)，發現程序對于出口高空泡份額的通道預測較準確。高功率棒束周圍的子通道普遍表現出高的出口空泡份額，包括邊通道和中間通道兩種。按照子通道周圍棒束功率的分布不同，又可將這兩種通道細分。表6列出選擇的具有代表性的子通道。

圖8 空泡份額預測值與實驗值的相對誤差

表7列出子通道79的不確定性分析結果。該結果表示，在考慮了表5所列的不確定性后，該程序對棒束實驗中子通道79出口空泡份額的計算，其結果為75.623%的可能性最大，有95%的可能性落在71.787%～79.966%之間，而這個判斷的可信度為95%。由圖2所示的功率分布不難看出，這5種子通道相比較，子通道17周圍的棒束功率最大，導致其空泡份額較高，符合本文分析的規律。

表6 8×8棒束實驗中典型高空泡子通道

表7 典型高空泡子通道的出口空泡份額計算不確定性范圍

計算了表6中5個子通道的不確定性，結果示于圖9。由圖9可見，高空泡子通道的出口空泡份額的計算不確定性差別很小，均在約-5.5%～6%之間。

圖9 典型高空泡子通道出口空泡份額的不確定性

子通道空泡份額的計算不確定性示于圖10。由圖10可見，不同子通道的計算不確定性差別很大。其中，邊角子通道的計算不確定性較小，約為±5.5%；水棒周圍不規則形狀的子通道的不確定性較大，約為±9%。因該類型子通道空泡份額的基礎小，受熱不均勻，且受熱較小，此時功率變化及交混量等因素的變化對其造成的影響不可忽略，因此不確定性較大。由此可見，該程序對高空泡子通道的出口空泡計算的不確定性較小，但由于該類子通道的空泡份額基礎較大，少許不確定性變動均有可能對堆芯安全造成威脅，因此應格外重視此類子通道空泡份額的不確定性的計算。

圖10 子通道空泡份額的計算不確定性

圖11 空泡份額的計算不確定帶與實驗值的比較

空泡份額的計算不確定帶與實驗值的比較示于圖11。該不確定帶滿足(95%，95%)的容忍限。然而，其結果并不能完全包絡實驗值，說明此次計算的誤差較大。其中，與實驗值差別較大的子通道有10、48、50、56和60，這幾個子通道的共同特征是其周圍棒束的加熱功率均不一致，即子通道周圍的受熱不對稱、不均勻。尤其是水棒周圍的子通道，其一面受到不斷加熱，而另一面又完全沒有受熱，這種極度的不對稱造成計算與實驗偏差較大。而子通道22、30、39、42的形狀規則，且是受熱對稱的通道，其計算的偏差較小且穩定。這說明在計算受熱不均勻的子通道時，子通道程序的計算能力還有待提高。

除此之外，造成不確定性范圍沒有包絡實驗值的原因還可能是：1) 基準工況的模型及模型選擇還存在不當之處；2) 模型參數的不確定性范圍無從得知，本文選定范圍并不合適，從而導致計算出的不確定性范圍與客觀事實不符。以上問題還需進一步探討。

4 小結

本文選擇BFBT的8×8棒束實驗作為計算對象，運用子通道分析程序對出口處的空泡份額進行了簡單的分析計算，針對子通道程序進行了不確定性分析，得到的結論如下。

1) 子通道程序對于邊角子通道或受熱不均勻子通道出口空泡份額的計算相對于實驗值的偏差較大。

2) 當邊界條件有微小變化時，出口空泡份額的變化量也較小。相比之下，入口焓的變化對結果影響最大。其次是出口壓力及熱流密度。對于出口空泡份額，子通道分析程序中的空泡漂移模型及修正種類對結果的影響較大，前者的最大影響高達25.4%，后者的最大影響為6.5%。在計算空泡份額時，這兩種模型的選擇非常關鍵。

3) 運用Wilks公式的原理確定能滿足容忍限的最少計算次數，且運用順序統計法對計算結果進行處理，得到計算值的不確定帶分布范圍。結果顯示，邊角子通道的計算不確定性較小，約為±5.5%；水棒周圍不規則形狀的子通道的不確定性較大，約為±9%。對于高空泡子通道的出口空泡份額，其不確定性在-5.5%～6% 之間，與其他子通道相比較小，但由于該類子通道的空泡份額基礎大，少許不確定性變動均有可能對堆芯安全造成威脅，因此應格外重視對該類子通道空泡份額的不確定性計算。所有子通道的不確定帶均未完全包絡實驗值，說明程序計算還存在較大偏差，需進一步改進。

參考文獻：

[1] Best estimate safety analysis for nuclear power plants: Uncertainty evaluation[R]. US: IAEA, 2008.

[2] BOYACK B E, CATTON I, DUFFEY R B, et al. Quantifying reactor safety margins, Part 1: An overview of the code scaling, applicability, and uncertainty evaluation methodology[J]. Nuclear Engineering and Design, 1990, 119(1): 1-15.

[3] NUPEC BWR full-size fine-mesh bundle test benchmark, Volume Ⅰ: Specifications[R]. US: NRC OECD Nuclear Energy Agency, 2005.

[4] 徐濟鋆. 沸騰傳熱和氣液兩相流[M]. 北京：原子能出版社，2001.

[5] WILKS S S. Determination of sample sizes for setting tolerance limits[J]. The Annals of Mathematical Statistics, 1941, 12(1): 91-96.

[6] GUBA A. Statistical aspects of best estimate method-Ⅰ[J]. Reliability Engineering and System Safety, 2003, 80(3): 217-232.