基于三維變分算法的蒸汽發生器模型同化方法研究

張　琪，楊淑涵，孫培偉，魏新宇，邱磊磊

基于三維變分算法的蒸汽發生器模型同化方法研究

張琪，楊淑涵，孫培偉*，魏新宇，邱磊磊

（西安交通大學核科學與技術學院，陜西 西安 710049）

蒸汽發生器在實際運行過程中系統參數會發生改變。為保證蒸汽發生器模型能夠真實反映蒸汽發生器的實際狀態，以蒸汽發生器水位、蒸汽壓力、一次側冷卻劑出口溫度以及循環流量等觀測變量，進行基于三維變分算法的蒸汽發生器模型同化多種蒸汽發生器觀測數據的方法研究。對蒸汽發生器仿真模型進行不確定性分析以確定參數調整范圍，基于SOBOL算法的敏感性分析結果選擇待優化參數組合，擬合關鍵參數與觀測變量之間的函數關系式作為觀測數據與模型的耦合方式，最后利用優化算法求解獲得了最優參數組合，實現了蒸汽發生器模型的數據同化。通過模型計算值與觀測數據對比，結果表明該方法可使蒸汽發生器的仿真精度明顯提高，證明了本文方法的正確性。

數據同化；蒸汽發生器模型；三維變分算法；不確定分析；敏感性分析

在壓水堆核電廠中，蒸汽發生器的安全運行和控制對整個核電廠的運行安全起著至關重要的作用。據統計，約25%的壓水堆緊急停堆均由不良的蒸汽發生器水位控制引起[1,2]。蒸汽發生器仿真模型是研究真實蒸汽發生器熱工水力動態特性的重要工具，廣泛應用于蒸汽發生器水位控制系統設計中，在水位控制系統設計中起到了巨大作用。因此在對發生器系統進行模擬時，精確的蒸汽發生器模型可以準確反映核電廠蒸汽發生器在實際運行時系統參數的動態變化過程，為操作員的調節提供指導和幫助。然而蒸汽發生器在長期的運行過程中系統參數會發生變化，因此不能忽視模型參數的不確定性對模型的仿真結果的影響，有必要研究數據同化方法，采用測量數據，自動校正模型的相關參數，保證模型能真實反映實際蒸汽發生器的實際運行狀態。這對設計合理的蒸汽發生器液位控制系統以及對于核電廠的安全穩定運行具有十分重要意義。

數據同化是一種將觀測數據和數學模型相結合，將觀測數據融入仿真模型對模型的初始條件和參數進行調整，以實現更好預測的方法。數據同化最初應用于氣象、水文等領域，近年來在核工程領域也有廣泛應用。Khalik[3]等建立了一個具有敏感性分析、不確定性量化和數據同化的集成框架，該框架可用于對具有大量輸入和輸出數據流的核模型進行最佳估計計算。Bouriquet和Argaud[4]等提出了一種對核質量進行最優評價方法，該方法對實驗數據和數值模型信息可對整合與評估后，進而對核質量進行調整，結果表明調整后的核質量冗余明顯降低。Eleveld[5]等為了改進核應急管理系統，對大氣彌散模型的數據同化方法進行了試點研究，結果表明核事故早期和晚期潛在污染區域的預測效果顯著提高。Bouriquet[6,7]利用三維變分算法和最佳無偏估計對測量網絡進行了優化設計，以重建900 MW壓水堆堆芯中子通量和功率場。Gong[8]等人提出了一種基于反距離擬合項的3D-Var數據同化方法，結果表明與已有的數據同化方法相比，該方法具有效果較好。Wan[9]等人提出了一種連續能量橫截面庫的核數據校正方法，通過核數據校正，實驗裝置的理論“模擬值”與相應的測量值吻合較好。

目前數據同化在核工程領域的成功應用主要集中于反應堆堆芯和核事故應急管理兩方面，在核動力裝置方向的研究很少。本文以蒸汽發生器為研究對象，基于三維變分算法，采用序列二次規劃算法對模型中的關鍵狀態進行優化調整，實現了針對蒸汽發生器的水位、循環流量、冷卻劑出口溫度和母管蒸汽壓力等多種觀測數據的同化，并驗證了蒸汽發生器模型同化方法的正確性。

1　模擬方法

1.1　蒸汽發生器模型

本文以“華龍一號”核電廠的蒸汽發生器為研究對象，由于無法獲得核電廠蒸汽發生器實際運行數據，因此依據真實蒸汽發生器的結構，采用系統程序RELAP5進行建模，并以此模型作為虛擬蒸汽發生器代替真實蒸汽發生器。

RELAP5是目前在輕水反應堆瞬態事故分析中最常用最典型的熱工水力分析程序，世界各國都對它做了很多驗證工作并生成了最少27篇報告。這些報告均使用真實有效的數據，根據實驗數據和RELAP5的計算結果比對，評估程序計算結果的精確性，使最佳估算程序在事故分析和安全分析中的應用得到普遍認可[10]。因此蒸汽發生器的RELAP5模型可以反映實際蒸汽發生器的特性，模型輸出可以代替真實的觀測數據。RELAP5模型節點如圖1所示，一次側冷卻劑進入入口水室后向上流動，流進U型管內，將熱量傳遞給二次側流體后，從出口水室流出；二次側給水從給水室向下流動，經過下降段，進入蒸汽發生器管束區被加熱而變成蒸汽，蒸汽繼續向上流動，在汽水分離器和干燥器中分離出飽和水的后從上封頭流出，經蒸汽母管到達汽輪機做功，被分離出的飽和水返回到再循環水室進行下一次循環流動。觀測數據選擇的是一次側冷卻劑出口溫度，二次側蒸汽發生器水位、循環流量和蒸汽母管壓力。

圖1　基于RELAP5的虛擬蒸汽發生器

1.2　蒸汽發生器仿真機

基于MATLAB/Simulink仿真平臺，依據質量、能量和動量守恒定理和可移動邊界理論，建立蒸汽發生器仿真模型，使用該模型實現與RELAP5模型產生的觀測數據進行同化方法研究。Simulink模型的節點如圖2所示，其中圖中節點編號P代表一次側，S代表二次側，M代表管壁。

圖2　Simulink蒸汽發生器模型

1.3　三維變分算法

三維變分算法是一種數據同化方法，它通過構建代價函數來描述模型計算值和觀測值的差異，對代價函數進行求解來獲得極小值。基于變分思想，把模型參數的優化問題轉化為一個極值求解的問題，在滿足約束的條件下，最小化模型計算值與觀測值之間的“距離”，使距離最小化的參數即為最優參數。三維變分算法代價函數如下：

式中：()——代價函數；

——待求解參數；

X——已知參數；

——背景誤差協方差矩陣；

——觀測數據；

——模型參數與觀測數據的映射關系；

——觀測誤差協方差矩陣。

1.4　研究流程

實際蒸汽發生器在運行過程中重點關注變量有水位、循環流量、冷卻劑出口溫度和母管蒸汽壓力等，選取RELAP5模型中這四個變量為觀測值對Simulink模型參數進行同化，在最小化代價函數之前，需要確定，，等變量，同化流程圖如圖3所示。

同化步驟為：

（1）模型的不確定性分析。假設模型參數的可調整范圍，采用廣義似然度不確定性分析，計算模型后驗分布和模型不確定度，確定模型參數調整的可行性；

圖3　同化流程圖

（2）敏感性分析。利用SOBOL算法對Simulink蒸汽發生器模型參數進行敏感性分析，為模型與觀測數據同化中的待優化參數的選擇提供依據。在敏感參數中，根據參數敏感度，實際物理意義以及盡量少參數的原則獲取參數組合；

（3）觀測算子的確定。使用Simulink模型，引入待優化參數擾動計算獲得觀測數據數值，通過函數擬合的方法獲得觀測算子；

（4）參數優化和效果檢驗。將計算獲得的協方差矩陣代價函數，采用最優化算法求解代價函數最小值。針對已確定的模型參數組合，利用Simulink模型計算獲得模擬結果與觀測數據進行對比和驗證。

2　模型不確定性分析

廣義似然度不確定性分析[11]（GLUE）是一種通用似然不確定性估計方法，該方法認為模型模擬結果的好壞應由模型參數組合決定，而不是單個參數。

2.1　模型參數選擇

模型計算中采用的換熱實驗公式都存在不確定度，常常可達±20%或±25%[12]，Simulink模型沸騰段采用的Chen換熱式的計算結果誤差甚至達到±50%[13]。阻力系數的不確定度對于模型計算結果的影響也不能忽視。根據蒸汽發生器的熱工水力過程以及蒸汽發生器Simulink模型的節點劃分，研究所選擇的主要模型參數、物理意義及其假設取值范圍如表1所示，由關系式計算所得數值乘以該修正因子就是模型參數的變化范圍。

表1　主要模型參數的意義及區間

2.2　不確定實驗設計

針對蒸汽發生器Simulink模型的換熱系數修正因子與阻力系數修正因子的不確定性研究，所采用實驗步驟如下：

（1）利用拉丁超立方（LHS）分層抽樣方法對每個模型參數進行均勻采樣，對各參數進行隨機組合，樣本組合為5 000組。

式中：——參數組對應的NS值；

——參數組的取值；

S——參數組在時刻對應的Simulink模型輸出值，m；

O——時刻RELAP5模型輸出值，m；

——仿真總步長，s。

（3）構建參數后驗分布和統計。排除步驟（2）中計算結果低于指定閾值的參數組，對高于指定閾值的參數組合進行歸一化，利用歸一化結果統計各參數的后驗概率密度分布。

（4）對步驟（3）中篩選出的參數組對應的模型模擬值進行排序，計算累計似然分布的5%和95%分位數，利用分位數結果獲得模擬水位的90%不確定性區間并計算不確定區間寬度，評估參數與仿真結果不確定性之間的關系：

式中：——不確定區間寬度，m；

——總時間步長，s；

——覆蓋度；

bracketed——位于90%置信區間的實測步長，s。

2.3　參數后驗分布與不確定區間

若某參數對納什系數無明顯影響，則參數的似然分布應與初始假設分布盡可能一致；若影響較大，則參數的似然分布與初始假設分布相差較大，呈正態分布。各模型參數的后驗分布如圖4（a）～圖4（h）所示。由圖可知，KR1、KR3、KR4和KR5的后驗分布呈近似均勻分布與所選的先驗分布基本相同，表明這4個參數對模型的輸出影響不明顯。KS3，KS4L，KS4R和KR2等參數后驗概率密度呈正態分布，且取值較為集中，說明這些參數對模型的水位輸出結果影響較大。Simulink模型模擬產生的水位90%不確定性區間如圖4（i）所示，由圖可知，蒸汽發生器Simulink模型模擬結果的90%不確定區間能夠完全包含（覆蓋率為100%，不確定性區間寬度為0.041 m）RELAP5模型的觀測數據，表明參數的假設取值范圍滿足要求。

圖4　模型參數后驗分布與不確定性區間

圖4　模型參數后驗分布與不確定性區間（續）

3　敏感性分析

3.1　目標函數

在進行敏感性分析前應根據研究模型自身特征選擇目標函數，選擇的目標函數應盡可能反映模型特征。選擇水位確定性系數，母管蒸汽壓力誤差系數，出口冷卻劑溫度誤差系數，循環流量誤差系數作為目標函數，四個目標函數可根據模型輸出與觀測數據的偏差的大小確定參數的敏感性。敏感性分析目標函數如表2所示。

表2　敏感性分析目標函數

續表

評價目標目標函數備注 蒸汽壓力 誤差系數蒸汽壓力相對誤差 溫度誤差 系數出口溫度相對誤差 循環流量 誤差系數循環流量相對誤差

注：為某時刻；obs() 為觀測水位；cal() 為模擬水位；obs() 為觀測蒸汽壓力；cal() 模擬蒸汽壓力；obs() 為觀測冷卻劑出口溫度；cal() 為模擬冷卻劑出口溫度；obs() 為觀測循環流量；cal() 模擬循環流量。

3.2　SOBOL指數法

式中：D——參數產生的方差；

D——參數和參數共同作用的方差；

D——參數參數和參數共同作用產生的方差；

1,2,3,…,m——個參數共同作用產生的方差。

將式（6）歸一化后可得到各參數和參數間共同作用的敏感性：

參數的各個指標：

一階敏感度：

二階敏感度：

總敏感度：

式中：S——參數單獨作用時的敏感度；

S——參數和參數相互作用的敏感度；

S——參數的主作用和參數與其他參數共同作用的敏感度之和；

——除了參數自身之外剩余參數的方差。

3.3　敏感性分析結果

以水位確定性系數為目標函數時，參數的一階敏感度與總敏感度如圖5（a）與圖5（b）所示，由圖知，敏感性參數為KS3、KS4L、KS4R和KR2，四個參數的一階敏感度分別為0.645、0.051、0.097以及0.165，總敏感度分別為0.689、0.072、0.114和0.188。以母管蒸汽壓力誤差系數為目標函數時，參數的一階敏感度與總敏感度如圖5（d）與圖5（e）所示，由圖知，敏感性參數為KS3，KS4L和KS4R。各參數的一階敏感度分別為0.311、0.099和0.231，總敏感度分別為0.431、0.321和0.445。以冷卻劑溫度誤差系數為目標函數時，各參數的一階敏感度與總敏感度如圖5（g）與圖5（h）所示，由圖知，敏感性參數為KS3，KS4L，KS4R。各參數的一階敏感度分別為0.425、0.163和0.285，總敏感度分別為0.424、0.233和0.354。以循環流量誤差系數為目標函數時，各參數的一階敏感度與總敏感度如圖5（j）與圖5（k）所示，由圖知，敏感性參數為KS3和KR2，各參數的一階敏感度分別為0.082和0.847，總敏感度分別為0.097和0.858。

由參數間的二階敏感度分析結果可知，各模型參數在圖中均有體現，表明兩個參數間的相互作用對蒸汽發生器模型的仿真結果影響十分顯著。

通過對四個目標函數下二階敏感度分析以及各參數的總敏感度可知，對于關鍵模型輸出影響明顯的為子KR2，KS3，KS4L，KS4R。四個敏感參數的組合作用對關鍵狀態量影響明顯，因此在同化實驗中需要對四個參數同時進行優化調整。

圖5　敏感性分析結果（續）

4　同化結果驗證與對比

4.1　蒸汽發生器數據同化過程

蒸汽發生器的數據同化由觀測數據、蒸汽發生器模型及其模擬結果的不確定性分析和觀測數據與模型的耦合方法等三個部分組成。關于觀測數據與模型的耦合需要通過一個觀測算子來建立聯系。這個觀測算子是模型參數關于可觀測狀態量的函數關系，應精確反映模型狀態，可以是簡單的線性函數或非線性模型。在蒸汽發生器模型模擬過程中，根據模型不確定分析結果，建立敏感性參數與關鍵狀態量線性函數關系。同化的實驗步驟如圖6所示。

4.2　穩態驗證

采用序列二次規劃算法（SQP）對代價函數在敏感參數的有效范圍內尋找最優解，在每次迭代過程中通過調整權重系數，來獲得使代價函數最小的模型參數組合。將每次調整后的參數作為新的模型初始參數帶入模型。5次同化后的模型穩態結果與穩態過程觀測數據對比如圖7所示。由圖可知，經過5次迭代后，水位由14.055 m變為13.847 m，實際觀測水位為13.843 m，相對誤差由1.568%降至0.002 8%。蒸汽壓力6.898 MPa升至6.955 1 MPa，實際觀測壓力為6.955 9 MPa，相對誤差由3.19%降低至0.012%。一次側冷卻劑出口溫度由289.49 ℃升至289.51 ℃，觀測出口冷卻劑溫度為289.61 ℃，相對誤差由0.041%降至0.035%。循環流量由2 817.07 kg/s升至3 198.58 kg/s，觀測循環流量為3 302.494 kg/s，相對誤差由15.24%降至3.42%。通過數據同化方法對模型參數進行調整，相對于同化前，水位、壓力、冷卻劑出口溫度和循環流量同化之后的精度分別提高了99.82%、99.62%、14.63%和77.55%。由此可見，穩態時，模型輸出的相對誤差明顯降低，同化前后的模型仿真精度提高明顯。

圖6　同化過程流程圖

4.3　瞬態驗證

在瞬態工況下，重點關注蒸汽發生器水位變化，利用瞬態過程中的水位觀測數據進行模型的回代檢驗。結合給水流量階躍下降2%和給水溫度階躍下降2 K兩種瞬態工況下RELAP5模型水位變化數據對同化后的模型進行外推檢驗，進一步驗證方法的有效性。回代檢驗和外推檢驗的結果如圖8所示。由圖8可知，對同化后的模型進行模擬檢驗，結果表明Simulink模型模擬水位與RELAP5模型觀測水位的決定系數R均在0.98以上，證明同化后的模型可以精確反映在瞬態工況時水位的變化趨勢。

5　結論

本文提出了一種含有不確定性分析、敏感性分析和數據同化的蒸汽發生器模型參數優化調整的方法。該方法實現壓水堆核電廠蒸汽發生器的觀測數據對蒸汽發生器仿真模型的關鍵參數進行調整，降低模型參數的不確定度，提高模型的仿真精度。

圖7　模型穩態結果對比

圖8　同化模型瞬態結果對比

本文得出的主要結論如下：

（1）對于蒸汽發生器，阻力系數和換熱系數是造成仿真不確定性的主要原因。

（2）該方法從模型自身出發，克服了依賴于專家經驗的手工法的局限性，利用觀測數據對與模型關鍵狀態量相關的敏感參數進行同化調整。

（3）對于瞬態工況，同化后的模型能精確反映水位變化趨勢，證明該方法可以大大提高模型的預測能力。

［1］ 張建民. 核反應堆控制［M］. 2版. 北京：中國原子能出版社，2016

［2］ 邱磊磊，張賢山，魏新宇，等. 自然循環蒸汽發生器的水位動態特性分析［J］. 核動力工程，2021，42（S2）：5-9.

［3］ Abdel-Khalik H，Turinsky P，Jessee M，et al. Stover T，Iqbal M.Uncertainty quantification，sensitivity analysis，and data assimilation for nuclear systems simulation［J］. Nuclear Data Sheets，2008，109（12）：2785-2790.

［4］ Bouriquet B，Argaud J P.Best linear unbiased estimation of the nuclear masses［J］. Annals of Nuclear Energy，2011，38（9）：1863-1866.

［5］ Eleveld H，Kok Y S，Twenhofel C J W.Data assimilation，sensitivity and uncertainty analyses in the Dutch nuclear emergency management system：a pilot study［J］. International Journal of Emergency Manage，2007，4（3）：551-563.

［6］ Bouriquet B，Argaud J P，Cugnart R.Optimal design of measurement network for neutronic activity field reconstruction by data assimilation［J］. Nuclear Instruments and Methods in Physics Research A，2012，664（1）：117-126.

［7］ Bouriquet B，Argaud J P，Erhard P，et al. Robustness of nuclear core activity reconstruction by data assimilation［J］. Journal of Power and Energy Systems，2012，6（2）：289-301.

［8］ Gong H l，Yu Y r，Li Q，Quan C y.An inverse-distance- based fitting term for 3D-Var data assimilation in nuclear core simulation［J］. Annals of Nuclear Energy，2020，141：107346.

［9］ Wan C H，Huang Y H，Zheng Y Q，et al. Nuclear-data adjustment based on the continue-energy cross-section library for the fast reactor［J］. Annals of Nuclear Energy，2020.124.

［10］蘇云，許以全，曹學武，等. 最佳估算程序在安全分析中的應用綜述［C］. 中國核學會&中國環境學會. 全國放射性流出物和環境監測與評價研討會論文集，2003：593-601.

［11］ Beven K，Binley A.The future of distributed models：Model calibration and uncertainty prediction［J］. Hydrological Processes，1992，6（3）：279-298.

［12］楊世銘，陶文銓. 傳熱學［M］. 4版. 北京：高等教育出版社，2006.

［13］于凱秋，孫立成，閻昌琪.“Chen”型沸騰傳熱計算方法用于小通道傳熱計算的適應性［J］核動力工程，2010，31（1）：24-27.

［14］羅鵬程. 武器裝備敏感性分析方法綜述［J］. 計算機工程與設計，2008：29（21）：5546-5549.

Research on Steam Generator Model Assimilation Method Based on 3D-VAR Algorithm

ZHANG Qi，YANG Shuhan，SUN Peiwei*，WEI Xinyu，QIU Leilei

（School of Nuclear Science and Technology，Xi’an Jiaotong University，Xi’an of Shaanxi Prov. 701149，China）

The system parameters of the steam generator change during actual operation. In order to ensure that the steam generator model can truly reflect the actual state of the steam generator, the method of assimilating the observation data of various steam generators based on the steam generator model based on the three-dimensional variational algorithm is carried out based on the observation variables such as the water level of the steam generator, the steam pressure, the primary side coolant outlet temperature and the circulating flow rate. The uncertainty analysis of the steam generator simulation model is carried out to determine the parameter adjustment range, the combination of parameters to be optimized is selected based on the analysis results of the SOBOL algorithm, the function relationship between the key parameters and the observed variables is fitted as the coupling mode between the observation data and the model, and finally the optimal parameter combination is obtained by using the optimization algorithm to solve, and the data assimilation of the steam generator model is realized. Comparing the calculated values of the model with the observational data, the results show that the simulation accuracy of the steam generator can be significantly improved, which proves the correctness of the proposed method.

Data assimilation; Steam generator model; 3D-Var algorithm; Uncertainty analysis; Sensitivity analysis

TL365

A

0258-0918（2023）05-1015-12

2022-09-09

自然科學基金面上項目（11875215）

張琪（1995—），男，陜西咸陽人，碩士研究生，現主要從事核反應堆控制與仿真方面研究

孫培偉，E-mail：sunpeiwei@xjtu.edu.cn