控制棒落棒時間分析軟件CRAC的研發與驗證

馬 超，周 毅，盧嘉川，蘆 韡，李 云，趙 波，陳 浩，馬 強，秦 勉，郭曉明

(中國核動力研究設計院 核反應堆系統設計技術重點實驗室，四川 成都 610213)

壓水堆緊急停堆時，控制棒組件依靠重力下落控制反應性，實現安全停堆。落棒時間計算是反應堆安全分析和燃料組件及相關組件設計的重要部分，關系到導向管及控制棒設計是否滿足落棒時間準則。針對落棒開展理論研究，研發落棒時間的分析工具，獲得落棒位移-速度-時間關系，對反應堆安全分析及燃料組件設計具有重要意義。

針對壓水堆落棒時間分析技術，國內外研究者開展大量研究工作，通過試驗與理論分析，認為落棒過程中控制棒驅動線受到的載荷主要為水力阻力和機械阻力。針對機械阻力，研究者利用有限差分或有限元法，求解驅動線各部件的振動方程，通過設置表面接觸條件，計算控制棒受到的非線性碰撞力與摩擦力[1-2]。水力阻力具體包括浮力、慣性力、局部形阻力和沿程阻力。針對水力阻力，研究者利用流體力學理論和經驗公式相結合的方法計算落棒過程中流體阻力[3]，還有研究者利用CFD軟件動網格技術，開展落棒過程中水力緩沖模擬[4-7]。上述研究工作考慮的落棒模型較為簡單，僅圍繞某一種載荷建模，如控制棒所受水力阻力或機械阻力，忽略或簡化考慮其他載荷因素，距離工程應用較遠。動網格技術需要耗費大量計算資源，只能構建簡單控制棒結構的仿真模型，如單根控制棒和水力緩沖部件，無法開展全尺寸控制棒驅動線模擬。

國內外核電企業采用理論和經驗關系式相結合的方式求解落棒過程中的水力阻力和機械阻力[8-10]，在試驗數據基礎上，形成針對特定堆型的落棒專用程序，開展落棒時間計算。本文介紹適用于華龍一號工程設計的落棒時間分析程序CRAC，滿足華龍一號落棒時間分析需求。

1 CRAC功能

CRAC通過合理簡化控制棒驅動線結構與載荷形式，建立動力學方程，基于流體力學壓降平衡關系、彈簧胡克定律及經驗機械阻力關系，求解落棒過程中的流體阻力與機械阻力，獲取反應堆正常工況及地震停堆工況，落棒總體時間及隨時間變化的控制棒位移、速度信息，適用于華龍一號控制棒驅動線的落棒時間計算。

2 理論模型

2.1 總體理論模型

控制棒驅動線結構如圖1所示，導向通道在軸向由下至上分為3部分：堆芯燃料組件部分、導向筒部件部分和控制棒驅動機構部分。反應堆正常運行狀態下，上部的驅動機構鉤爪提升驅動桿，將控制棒組件大部分停留在堆芯外部，控制棒棒束保持在導向筒部件中；需要停堆時，鉤爪部件斷電，松開驅動桿，控制棒驅動線依靠重力下落到燃料組件中，實現停堆。驅動桿、控制棒為細長結構(長徑比大于100∶1)，控制棒驅動線與通道之間間隙較小，限制了徑向運動，因此在分析落棒過程中，近似為一維軸向運動。

圖1 反應堆控制棒驅動線結構Fig.1 Structure of reactor control rod drive line

圖2示出控制棒驅動線受力示意圖，針對控制棒軸向運動建模，所有作用力都作用在重心，建立運動方程如下：

圖2 控制棒驅動線受力示意圖Fig.2 Force analysis schematic of control rod drive line

(1)

2.2 流體阻力計算模型

控制棒在驅動機構、導向筒及燃料組件3部分通道中運動。由于在導向筒中流場較為寬闊，水力阻力較小，故忽略導向筒內作用在行星架與控制棒上的水力阻力。運動部件所受流體阻力，主要考慮控制棒在燃料組件中及驅動桿在驅動機構中受到的流體阻力。引入如下假設：1) 流體不可壓縮；2) 落棒過程中在燃料組件導向管和緩沖段內的水溫為常數；3) 在每個計算時間步內，固定壁面摩擦系數與運動壁面摩擦系數為定值(兩者選取不同參考速度)。

2.2.1控制棒在燃料組件內所受到的水力阻力 圖3示出控制棒插入到導向管和緩沖段示意圖，采用集總參數法描述流道速度信息，如表1所列。由于堆芯流量可能在較大范圍內變化及控制棒組件機械緩沖回彈現象存在，各流道速度方向可能發生變化。在水力阻力建模前，需定義速度標志位，保證對流體動力學方程求解的正確。

圖3 控制棒插入到導向管和緩沖段示意圖Fig.3 Schematic of control rod dropin thimble guide and dashpot

表1 各通道流體速度定義Table 1 Definition of fluid velocity direction in different flow channels

1) 控制棒在導向管內運動

在燃料組件中所受水力阻力，由壓差引起的作用力和沿程阻力組成：

FFA=NCR(SCR(pB-p2)+lABCRτCR)

(2)

式中：FFA為控制棒組件在燃料組件內受到的水力阻力；NCR為控制棒組件中棒總數；SCR為單根控制棒橫截面積；pB為控制棒下端界面處壓強；p2為導向管頂端外部壓強；lA為控制棒插入導向管軸向長度；BCR為單根控制棒濕潤周長；τCR為控制棒在導向管內受到的平均流體切應力。

連續性方程為：

vSCR-vASA-4vDHSDH=0

(3)

4vDHSDH=vTGSTG

(4)

式中：STG為導向管圓形通道過流面積；SA為導向管環形通道過流面積；SDH為單個流水孔過流面積(周向分布共4個)。流水孔外側和控制棒下端之間的流體壓差為：

(5)

式中：p1為流水孔外側的流體壓強；ξDH為流水孔局部阻力系數；λTG為導向管圓形通道的沿程阻力系數；DTG為導向管通道水力直徑；lTG為導向管總長度。

在導向管環形通道內，控制棒下端與導向管上端之間壓差表示如下：

(6)

λ按Colebrook公式計算[11]：

(7)

式中：Ks為壁面絕對粗糙度；Re為雷諾數，Re=vrefDρ/η，vref為參考速度，η為流體黏度；D為水力直徑。

從流水孔外側到導向管頂端壓降Δp12=p1-p2，由燃料組件特性、進入燃料組件流量及流體溫度共同確定，因此Δp12視為已知量，根據壓降平衡關系：

Δp12=(p1-pB)+(pB-p2)

(8)

在確定每一時間步v(t)條件下，聯立式(3)～(6)及(8)可確定vA(t)，從而確定當前時刻控制棒壓差(pB-p2)和切向應力分量：

(9)

根據式(2)獲得控制棒受到的水力阻力FFA。

2) 控制棒在緩沖段內運動

當控制棒插入緩沖段時，控制棒組件受到的水力阻力FFA由式(10)計算，此時導向管內環形流道長度為lA=lTG。

(10)

式中：lB為控制棒插入緩沖段內長度；CCRB為緩沖段環形流道中控制棒壁面摩擦系數。此時，連續性方程如下：

4vDHSDH+vASA-vBSB=0

(11)

vSCR+vDPSDP+vBSB=0

(12)

式中：SB為緩沖段環形流道過流面積；SDP為緩沖段圓形流道過流面積。控制棒下端與導向管頂端之間壓降分為兩段：

(13)

(14)

式中：pDH為流水孔內側流體壓強；CDP為環形流道緩沖段壁面側摩擦系數；ξ1BS為環形流道截面突然擴張的局部壓力損失系數，即從SB突擴至SA的壓力損失系數。由導向管緩沖段對應的外部流道壓降平衡關系為：

pB-pDH=(pB-p0)+

(p0-p1)+(p1-pDH)

(15)

式中，p0為軸肩螺釘阻尼孔底端壓力。

(16)

(17)

式中：lDP為緩沖段總長；λDP為緩沖段圓形通道沿程損失系數；λSH為軸肩螺釘阻尼孔中沿程損失系數；DDP為緩沖段圓形通道水力直徑；DSH為軸肩螺釘孔水力直徑；lSH為軸肩螺釘孔深；ξSH為軸肩螺釘阻尼孔局部壓力損失系數；SSH為軸肩螺釘孔過流面積。因此式(15)可轉化為：

(18)

式中，Δp01=p0-p1視為已知。

在已知v(t)情況下，聯立式(13)、(14)、(16)和(18)可確定vA(t)、vB(t)，進而確定壓差(pTG-pA2)，通過式(10)計算控制棒插入緩沖段時受到的水力阻力FFA。

2.2.2驅動桿在驅動機構內受到的水力阻力 驅動桿在驅動機構中受到的水力阻力建模過程與2.2.1節內容相似，通過求解連續性方程和動力學方程，可得驅動桿受到的水力阻力FCRDM。

2.3 機械阻力計算模型

在控制棒驅動線落棒過程中，驅動桿與控制棒會發生變形，與導向通道發生隨機碰撞接觸，機械運動行為復雜，直接針對驅動線相關設備結構形式建模，計算模型將非常龐大且求解困難。因此在CRAC研發過程中，考慮基于控制棒驅動線相關試驗結果與工程經驗，建立機械阻力相對于錯對中量及插入深度變化的計算模型：

Ff=F(d,z)+F0

(19)

式中：函數F(d，z)輸入變量包括上部導向組件與隔熱套之間錯對中量d及落棒位移z，F(d，z)是在參考試驗數據條件下，考慮工程保守性建立的插值函數；F0為對中條件下的機械摩擦力。

3 CRAC計算流程

圖4示出CRAC計算流程，通過輸入整個系統初始參數，包括結構尺寸參數、冷卻劑物性參數、控制棒初始位置與速度、各通道中流體初始速度和壓強，確定控制棒驅動線初始加速度。將時間步長設為Δt，第i時間步對應時刻ti=ti-1+Δt，該時刻所有載荷與流體狀態將根據時刻ti-1的特性確定。每一時間步求出控制棒驅動線所受合力，求得加速度及其他下落特性參數。重復上述過程，直至控制棒組件完全插入燃料組件，實現對落棒過程的模擬計算。

圖4 CRAC計算流程圖Fig.4 Calculation flow chart of CRAC

4 CRAC驗證

4.1 落棒試驗簡介

CF2燃料組件是由中國核動力研究設計院自主設計研發，是第三代核電華龍一號核心部件。為滿足第三代核電燃料組件剛度要求，CF2在導向管緩沖結構方面采用管中管新型設計。該燃料組件將應用于華龍一號國外首堆巴基斯坦卡拉奇2號機組的首循環。本次落棒試驗中回路參數為：壓力15.5 MPa，溫度315 ℃，工況列于表2[12]。試驗測量得到：控制棒下落至緩沖段入口的落棒時間(T5)；控制棒下落至全插入位置的落棒時間(T5+T6)；控制棒下落速度、位移與時間關系。落棒時間測量誤差為±0.01 ms，落棒速度測量系統相對誤差為±0.28%。

4.2 計算結果與試驗值對比

根據試驗條件及控制棒驅動線結構參數，利用CRAC進行表2所列5種工況下的落棒模擬計算。圖5示出G1和G5工況下速度-時間和位移-時間曲線。由圖5可知：CRAC計算結果與試驗結果變化趨勢一致，但計算最大速度值要比試驗值略低；CRAC計算的位移曲線趨勢較為平緩，落棒時間結果要比試驗值大。因為CRAC中機械摩擦阻力及水力局部壓降系數均為經驗函數或經驗系數，充分考慮落棒時間計算保守性，保證軟件計算值大于真實值。

圖5 工況G1(a)、G5(b)條件下的落棒結果對比Fig.5 Comparison of control rod drop result for G1 (a), G5 (b) conditions

表2 落棒試驗工況Table 2 Condition of drop experiment

落棒試驗獲得的速度曲線在控制棒進入緩沖段之前達到最大，此時曲線出現明顯的突然上升、然后下降過程。產生該現象的可能原因是，控制棒下落的過程中，子彈頭部位會產生局部過壓，在導向管中運動過程中，距離流水孔較遠時，過壓產生的水力阻力壓差增大較為平緩連續，而在控制棒子彈頭接近流水孔時，流水孔處內外壓差較大，流水孔起到泄壓作用，使得控制棒子彈頭部位過壓在流水孔附近產生先降低(泄壓)再升高(進入緩沖段)的變化，因此速度曲線在此處先迅速升高再迅速下降。由于CRAC采用一維集總參數模型，將控制棒簡化為圓柱體，無法考慮子彈頭部位的壓力場與速度場急劇變化現象，因此獲得的速度-時間曲線是逐漸變化，沒有突然上升再下降的過程。

表3列出試驗數據(T5，T5+T6)與計算結果對比，計算值均大于試驗值，實現對試驗數據的包絡，滿足反應堆安全分析保守性要求。落棒時間計算值同試驗值符合性較好，最大相對誤差為15.45%，證明CRAC對華龍一號堆型落棒時間計算的適用性以及準確性。

表3 落棒時間計算與試驗數據對比Table 3 Comparison of calculated and measured values for control rod drop time

5 結論

本文介紹了CRAC的功能、理論模型、計算流程，并利用試驗數據對CRAC的計算準確性進行驗證，計算結果與試驗數據相對誤差在15.45%以內，兩者符合較好，且計算結果均大于試驗結果，表明軟件計算精度與保守性能滿足華龍一號堆型安全停堆時間準則分析需求。