典型壓水堆小破口LOCA參數重要度排序表開發

熊青文，黃 濤，茍軍利，杜 鵬，鄧 堅，袁 鵬，周佳樾，胡文楨

(1.中國核動力研究設計院 核反應堆系統設計技術重點實驗室，四川 成都 610213；2.西安交通大學 核科學與技術學院，陜西 西安 710049)

最佳估算加不確定性(BEPU)分析方法是國際原子能機構(IAEA)推薦用于核電廠安全分析和執照申請的先進方法[1]，該類方法的一般步驟為：首先基于確定的目標輸出(FOM)建立現象識別排序表(PIRT)識別重要現象、模型或參數，并以帶有真實不確定性分布的輸入表示這些不確定性源；而后使用最佳估算程序建立電廠的模型，并執行多次程序計算將輸入參數的不確定性傳播至FOM；最后基于程序計算結果量化FOM的不確定性，并執行相應的不確定性分析和敏感性分析。由于核電廠的結構十分復雜，因此在模擬電廠時會涉及數量龐大的輸入參數，如初始/邊界條件、材料物性、狀態參數、本構模型等[2]。在BEPU分析中考慮所有參數的不確定性是不現實的，因此需要基于一定的方法篩選和識別在分析工況中比較重要的參數以簡化分析。其中，PIRT是最常使用的一種方法。通常而言，PIRT的建立需要使用一定的方法，如AHP[3]或MCDA[4]等，并在專家經驗的指導下使用定性的方法表征現象或參數的重要度，如使用H、M和L分別表征高、中、低重要度[5]。PIRT能初步識別事故分析過程中的重要現象、模型或參數，但是在實際應用中會存在一些問題。首先，PIRT建立的成本很大。由于不同反應堆在不同工況下的事故序列和熱工水力現象存在區別，因此理論上針對不同的反應堆堆型和不同的事故工況均要建立對應的PIRT，導致很多特定的反應堆不存在完善的PIRT，難以開展BEPU分析而使用保守參數進行保守分析。此外，PIRT的建立高度依賴于專家經驗，在專家經驗不足的情況下往往采取的措施是參考類似反應堆的已有PIRT，這種方式存在可能將不重要的參數考慮在內，而將重要參數遺漏的問題。

核電廠安全分析中失水事故是最為關注的工況之一，相比于大破口失水事故而言，小破口失水事故(SBLOCA)的BEPU分析研究較少，相關的PIRT也存在不多。為了能在缺乏PIRT的情況下執行準確的參數重要度排序計算，本文以典型三回路壓水堆(PWR)不同破口尺寸下的SBLOCA為對象開展相關工作，選定的FOM為堆芯活性區最低坍塌液位(MCL)。首先，分別在冷管段橫截面積的0.2%、0.4%、0.6%、0.8%和1.0%的5個破口尺寸下開展SBLOCA的瞬態基準計算，確定工況中的階段劃分，并確定堆芯MCL會出現于噴放和環路水封存在及清除階段兩個階段。隨后，使用優化矩獨立全局敏感性分析方法對5個破口工況中兩個階段分別執行定量敏感性分析，并將參數的重要度排序轉換為Savage分數，按照Savage分數將所有輸入參數進行重要度分組，從而得到PWR SBLOCA的參數重要度排序表。

1 優化矩獨立方法介紹

本文使用的敏感性分析方法為矩獨立全局敏感性分析方法[6]，其旨在評估輸入參數對目標輸出概率密度函數(PDF)的影響，通過固定輸入參數的取值并計算目標輸出的無條件和條件PDF間的偏移量來定量評估參數的敏感性度量。假設某個模型Y=g(X)存在k個輸入參數，第i個輸入參數的矩獨立敏感性度量δi可表示為：

(1)

(2)

式中：fXi(xi)為第i個輸入參數的概率分布，xi為第i個輸入參數Xi在計算中的取值；s(Xi)為固定第i個輸入參數情況下輸出PDF的偏移量；fY(y)、fY|Xi(y)分別為Y的無條件和條件PDF，y為輸出參數Y在計算中的取值。

由于計算參數的矩獨立敏感性定量涉及到兩個嵌套積分計算，因此總體的計算量十分巨大。考慮到核電廠SBLOCA工況的模擬特別耗時，一次計算需耗費數小時，因此難以采用常規抽樣方法直接計算各輸入參數的矩獨立敏感性度量。為了能使計算成本特別小且又可相對精確地計算輸入參數的敏感性度量，本文使用一種優化的矩獨立全局敏感性分析方法。首先，為了降低積分計算的計算量，使用五點高斯求積方案替代積分計算[7]。基于高斯求積方案，可將矩獨立敏感性度量表示為：

(3)

式中：ωij為第i個輸入參數按照其分布類型確定的5個高斯權重值中的第j個值；Xij為第i個輸入參數的第j個高斯點取值。因此求解度量δi的關鍵在于求解s(Xij)，即輸出的無條件PDF與將第i個輸入參數固定于其第j個高斯點時輸出條件PDF間的偏移量。

Liu等[8]提出了一種基于輸出Y的條件和無條件累計分布函數(CDF)計算s(Xij)的方法，計算得到輸出的CDF即可快速計算s(Xij)。根據CDF的定義可表示為：

FY(y)=P{Y≤y}=P{g(X)≤y}=

P{g(X)-y≤0}=P{z(X,y)≤0}=

Pf{z(X,y)}

(4)

式中：z(X,y)=g(X)-y為定義的新函數；Pf為失效概率。因此，可將求解模型Y=g(X)的CDF轉換為求解模型z(X，y)的失效概率。優化方法中，使用四階矩估計方法[9]和Pearson系統[10]來求解模型的失效概率。

由于采用五點高斯求積方案簡化計算，因此每個參數僅需計算5個高斯點處對應的函數輸出值，同時執行1次所有參數的名義值計算，因此該方法需要調用的模型計算次數最多為5k+1次。傳統計算矩獨立敏感性度量所需的計算次數常為輸入參數數量的數萬至數百萬倍，相比于傳統基于隨機抽樣的矩獨立敏感性度量計算方法而言，本文提出的方法能夠以極小的計算量達到類似的計算精度。使用失水事故試驗(LOFT)中的LP-02-6工況對提出的方法進行驗證[11]，驗證計算中使用了15個重要輸入參數，分別使用提出方法和傳統抽樣方法計算輸入對包殼峰值溫度的影響，其中抽樣方法抽樣計算次數為1.001×106，而優化方法總計算次數僅為81。兩種方法的計算結果如圖1所示，結果表明，在計算量縮小近1萬倍的情況下，對于大部分輸入參數該提出方法與傳統方法計算結果的相對偏差在0.8%～5.2%內。

圖1 LOFT LP-02-6工況敏感性分析結果

2 SBLOCA瞬態分析

本文考慮的工況為典型PWR 5個不同破口尺寸下的SBLOCA，首先執行所有參數取名義值的基準工況計算，并基于計算結果分析SBLOCA瞬態。計算結果表明，SBLOCA的事故進程可根據現象分為4個階段，這4個階段及其主要現象介紹如下。

1)噴放階段

小破口出現后，冷卻劑快速從破口處流失，導致一回路壓力快速下降。但由于反應堆此時還未觸發停堆狀態，主泵還在持續運轉，一回路的冷卻劑處于受迫循環狀態。當反應堆停堆后，保守假設主泵立即停止運行，主泵在惰轉一段時間后停止，此時一回路內冷卻劑進入自然循環階段。

2)自然循環階段

主泵停止轉動后，一回路冷卻劑處于自然循環階段。此時堆內冷卻劑依靠密度差進行流動，自然循環階段持續時間往往較短。

3)環路水封存在及清除階段

當一回路壓力降至飽和壓力時，一回路內冷卻劑開始沸騰，產生大量蒸汽。由于破口尺寸較小，蒸汽不能及時從破口排出，導致積聚在上腔室和蒸汽發生器(SG)U型管內，蒸汽的積聚又阻礙了冷卻劑的流動，進一步導致換熱惡化，這個現象被稱為環路水封。環路水封的出現會導致堆芯液位快速下降。隨著一回路冷卻劑從破口的不斷流出，環路水封被逐漸清除，一回路壓力恢復下降，且堆芯液位開始被安注系統恢復。

4)長期冷卻階段

環路水封清除后，堆芯液位開始恢復，直至完全淹沒，過冷冷卻劑能持續穩定地帶走堆內的衰變熱。

對于破口尺寸較小的SBLOCA工況(如0.2%破口工況)而言，其堆芯活性區MCL出現在噴放階段，主要原因是閃蒸現象導致的堆芯局部沸騰及破口的冷卻劑流出。而對于破口尺寸較大的SBLOCA工況而言，其堆芯活性區MCL出現在環路水封存在及清除階段，主要原因是堆內傳熱惡化導致的堆芯強烈沸騰及破口的冷卻劑流出。由于兩種情況下堆芯活性區MCL出現的原理不同，后續敏感性分析應分為兩個階段獨立開展，以評估各輸入參數在噴放階段和環路水封存在及清除階段對堆芯活性區MCL的影響。

3 敏感性分析

由于缺乏PWR SBLOCA相關PIRT，因此擬參考類似研究中建立的PIRT或使用的重要參數，并結合目標電廠補充相關參數。諸多針對SBLOCA的PIRT中，西屋公司針對AP600開發的PIRT包含了典型能動式壓水堆和非能動式壓水堆中的主要現象，具有較高的可信度，因此本文主要參考的PIRT為AP600核電廠的SBLOCA PIRT[12]，但不包含非能動設備相關的現象或參數。此外，還參考了文獻[13-15]中對SBLOCA工況開展BEPU分析時使用的重要參數。基于參考類似PIRT共識別了34個輸入參數，此外，通過分析基準工況計算中出現的重要現象，補充了可能存在重要影響的模型或參數共20個輸入，因此在敏感性分析計算中共考慮了54個輸入參數。根據提出的優化矩獨立全局敏感性分析方法，54個輸入參數執行全局敏感性分析需進行217次程序計算。因此，使用所有54個不確定性輸入參數執行全局敏感性分析，于5個破口尺寸工況各執行217次程序計算，并進一步計算54個輸入參數在噴放階段和環路水封存在及清除階段對堆芯活性區MCL的矩獨立敏感性度量。計算結果如圖2所示，由于0.2%破口尺寸下沒有環路水封存在及清除階段，因此沒有相應的敏感性分析結果。

圖2 0.2%、0.4%、0.6%、0.8%、1.0%破口尺寸工況敏感性分析結果

根據Pareto法則[16]，復雜模型可被小部分重要的參數和大部分影響較小的參數表征。圖2示出了不同破口尺寸和不同事故階段中各輸入參數對堆芯活性區MCL的影響，表明大部分影響集中于少數重要輸入參數，該結果符合Pareto法則。考慮到缺乏SBLOCA的PIRT，因此可基于敏感性分析結果建立參數重要度排序分組表，可為PWR SBLOCA工況BEPU分析提供參考和指導。需要說明的是，由于大部分影響集中于少數重要輸入參數，若按照敏感性度量進行參數排序和重要度分組，會出現重要度排序為高的參數數量較少的情況。而且不同破口尺寸工況計算得到的各參數的敏感性度量無法直接比較，因此可將參數的重要度排序轉換為Savage分數，然后將所有參數進行重要度定性均勻分組，分別為高(H)、中(M)、低(L)和無影響(N)。Savage分數的定義為：

(5)

式中：N為所有輸入參數的數目；Ri為第i個輸入參數的秩排序。

根據計算結果，噴放階段有16個參數沒有影響，其重要度排序為54，分組為N。剩余38個參數進行均勻分組，即排序1～13的參數分組為H，排序14～26的參數分組為M，排序27～38的參數分組為L。同理，在環路水封存在及清除階段有7個參數對MCL沒有影響，其分組為N，并對剩余47個參數進行均勻分組。根據最終對結果的處理，可得到PWR在SBLOCA工況下的參數重要度排序表(表1)。由表1可知，噴放階段對MCL影響較大的參數主要為初始功率、初始一回路壓力、臨界流模型以及一回路中的流動阻力。而在環路水封存在及清除階段，除破口冷卻劑損失外，堆芯液位下降還與堆內沸騰有關，因此界面阻力模型、各種換熱模型以及輔助給水系統延遲也對MCL有較大的影響。對比發現，提出方法計算得到的敏感性分析結果與參考PIRT較為相似，此外亦與現象機理具有高度一致性，因此可說明該敏感性分析結果的可信度。需要補充說明的是，該計算得到的PIRT能適用于典型能動式PWR SBLOCA工況的現象識別，亦可用作非能動式PWR SBLOCA工況的現象識別參考。

表1 PWR SBLOCA工況參數重要度排序分組表

4 結論

本文研究了在缺乏相應PIRT的情況下確定各輸入參數的重要度，選取的研究工況是典型三回路PWR在冷管段發生破口尺寸分別為0.2%、0.4%、0.6%、0.8%和1.0%下的SBLOCA。首先基于基準計算結果，結合已有的SBLOCA PIRT，篩選了可能對FOM，即堆芯活性區MCL，具有影響的54個不確定性輸入參數。然后，使用一種優化矩獨立全局敏感性分析方法，計算得到了不同破口尺寸下各輸入參數對FOM的敏感性度量和重要度排序。將參數的重要度排序轉換為Savage分數，按照Savage分數定性地將所有輸入參數進行重要度分組，從而得到了PWR SBLOCA的參數重要度排序表。本文提供了一種在缺乏PIRT時識別重要參數的方法。計算結果表明，該文使用的優化矩獨立方法能以較小的計算成本完成全局敏感性分析，并能準確地識別對FOM沒有影響的輸入參數。基于矩獨立敏感性分析的計算結果，建立了PWR在SBLOCA下的參數重要度排序表，該表能夠在SBLOCA相關分析缺乏PIRT的情況下提供一定的指導和參考價值。敏感性分析結果表明，噴放階段對MCL影響較大的參數主要為初始功率、初始一回路壓力、臨界流模型以及一回路中的流動阻力。而在環路水封存在及清除階段，除破口冷卻劑損失外，堆芯液位下降還與堆內沸騰有關，因此界面阻力模型、各種換熱模型以及輔助給水系統延遲也對MCL有較大的影響。