控制棒下落时间计算模型

刘 新 陈先龙 张 修 于晓雷

（深圳中广核工程设计有限公司 深圳 518172）

控制棒下落时间计算模型

刘 新 陈先龙 张 修 于晓雷

（深圳中广核工程设计有限公司 深圳 518172）

反应堆驱动线是反应堆控制的“生命线”，执行重要的安全功能。控制棒下落时间不仅是驱动线设计考核的重要指标，也是核电站安全分析的关键参数。基于压水堆驱动线的结构特性，对控制棒组件在下落过程进行受力分析，提出落棒时间计算的数学模型，采用有限差分法求解运动方程，并且在此基础上编制计算程序。计算结果最终表明：程序计算值与岭澳核电站二期机组热态落棒试验值符合较好，基于此计算方法编制的落棒时间计算程序可以用于压水堆停堆时的安全评价。

控制棒组件，驱动线，落棒时间

控制棒组件在反应堆中的落棒时间对于反应堆的安全运行至关重要。反应堆驱动线由运动组件和固定组件两部分组成。运动组件包括控制棒棒束组件(Rod Cluster Control Assembly, RCCA)和驱动杆，固定组件是运动组件的下降通道，主要包括：控制棒驱动机构(Control Rod Drive Mechanism, CRDM)、热套管、控制棒导向筒组件和燃料组件导向管。反应堆驱动线各部件结构复杂，落棒时间计算需要考虑反应堆内冷却剂的流速、温度和压力的影响，也要考虑外界因素（如地震）的影响。

控制棒在下落过程中，由于通道流体流速的变化，导致受到的水力学力（水力摩擦力和水力压力）会不断的变化，而且变化值与下落位置是非线性关系，从而增加了计算落棒时间的难度[1]。窦一康等[2]提出了地震工况下用有限差分法对振动方程作离散，对时间历程则用Wilson-θ法进行逐步分析，对由碰撞引起的非线性采用修正的Newton法作迭代处理的计算方法。但是，常微分方程的 Wilson-θ求解方法比较复杂，因此，我们采用有限差分方法，将下落时间进行离散，同时利用圆管中流体压差计算的原理，避开了复杂的常微分方程计算[3]。此外，目前国内外有关控制棒下落时间计算的软件大多采用Fortran面向过程的编程语言，不利于代码的阅读以及程序的扩展，本文在总结出更精确、更简单计算模型的基础上，采用面向对象编程语言C++编写出了控制棒落棒时间计算程序[4]。

1 理论模型

控制棒运动组件不仅受到重力还受到阿基米德浮力、水力压力、水力摩擦力等，其中有些力随时间进行变化，为了计算落棒时间，我们需要对计算模型进行简化假设，并且要把运动过程按照时间步长进行离散，假设在单个时间步长过程中力是恒定的，然后对时间进行积分计算。

1.1 简化假设

(1) 在控制棒导向筒内，作用在星型架和控制棒的水力学力被忽略；

(2) 导向管内的水温为常数且等于堆芯中冷却剂的平均温度；

(3) 在计算的任何时刻，流体应为单向的。

1.2 控制棒受力分析

控制棒下落过程受到的力包括恒力和变力（如图1所示），其中，恒力有重力F1、浮力F2、额定机械摩擦力F3、控制棒导向筒在连续导向段的水力-机械摩擦力F4，变力包括控制棒在导向管中受到的水力学力F5、驱动杆在驱动机构中受到的水力学力F6、地震工况下的附加力F7。

图1 控制棒下落过程受力分析Fig.1 Force analysis of control rod in dropping.

根据牛顿第二定律，有如下动力学方程：

式中，M为控制棒运动组件的质量，kg；a(t)为t时刻控制棒运动组件的加速度，m·s-2。

由式(1)可以看出，需要计算下落时间的关键因素就是变力的计算，因为控制棒运动组件是在流体中运动，我们可以根据流体在圆管中的流动模型以及圆柱体插入圆管中的流动模型分别进行分析，详细模型说明见§1.3和§1.4。

1.3 流体在圆管中流动时的模型

假定流体不可压缩，不考虑重力的影响，流体沿管道流动，如图2所示。流动截面上的静压和壁面摩擦力均匀分布，则动量方程为：

式中， Cf为壁面摩擦系数；U为流体的流速，m·s-1；Pe为圆管的湿周，m；S为管道的截面积，m2。

图2 流体在圆管中流动模型Fig.2 Fluid model of flows in a pipe.

则动量方程[5]可写成：

进行空间积分：

1.4 圆柱体在圆管中运动时的流动模型

1.4.1 压降计算

圆柱体在管道中运动，如图3所示。假设圆柱体在管道中的下落速度为V，根据式(2)得出流体此时的动量方程：

式中，Cf′为圆柱体的壁面摩擦系数；Pe'为圆柱体的湿周，m；Sa为圆柱体与壁面环腔之间的截面积，m2。

通过对动量方程积分得：

图3 圆柱体在圆管中流动模型Fig.3 Fluid model of flow in a pipe with a mobile cylinder.

1.4.2 圆柱体受力分析

流体和圆柱体间的水力学摩擦力通过壁面摩擦系数Cf′来确定：

由于流体流动产生的压降而作用在圆柱体顶部和底部的压差阻力根据式(8)进行计算：式中，SCR为圆柱体的截面积，m2；zΔ为控制体的步长，m。

1.5 控制棒组件下落过程中变力分析

控制棒组件在下落通道中受到的流体阻力计算主要的难点为驱动机构中驱动杆受到的流体阻力和控制棒在燃料组件导向管中受到的流体阻力的计算，根据§1.3和§1.4的物理模型，我们采取下面的思路进行计算。

1.5.1 驱动机构区域

重点在于计算驱动杆在驱动机构区域流体流速的计算，把驱动结构分为三段求解，分别为行程套管区、钩爪区、热套管区。在三个作用区域分别利用流体的质量守恒方程和§1.4给出的压差和受力计算方程联立得出非线性方程组，然后采用牛顿迭代法进行求解，最后计算出驱动杆在行程套管区、钩爪区、热套管区三个区域受到的流体水力作用力。

1.5.2 燃料组件导向管区域

重点在于计算燃料组件导向管中控制棒下方流体的流速以及控制棒与导向管环形通道中流体流速。根据§1.3和§1.4压差的计算方程，得出导向管内流体的压差的方程，然后根据燃料组件外部的压差与内部压差的关系，计算出导向管内各处流体的流速，最后计算出控制棒在导向管内受到的流体水力作用力。

2 算法分析

从初始时刻t(0)开始计算，按时间步长迭代计算。初始状态工况的描述如下：

(1) 移动组件的初速度为零；

(2) 导向管内的流体流量是恒定的；

(3) 在管座/热套管/驱动机构部位的流体流量为零；

(4) 控制棒被提出或被部分提出。

在 t(i) 时刻，我们计算每个产生作用的力的值。我们从中得出加速度，并且假设在单个时间步长中力是恒定的，得出 t(i+1) 时刻的速度和行程。基于这些值，重新计算变力。

由上可以看出，在每个时间步长之后，需要重新计算的力有导向管内控制棒的水力学力F5、驱动机构部位驱动杆的水力学力F6、附加力F7（地震工况下才需要计算）。我们计算控制棒从任意给定位置依靠重力下落到导向管缓冲段入口的时间以及星型架底端弹簧刚好接触上管座上表面的时间。

在地震工况下，需借助试验数据通过插值的方法得出t时刻控制棒组件受到的地震附加力[6]，但尚未进行相关的模拟振动试验，所以暂时预留接口。

利用面向对象的思想，对控制棒下落过程进行模型抽象化，创建阻力模型管理器，把下落过程分解出的阻力计算模型逐个添加到阻力模型管理器中。计算流程如图4所示。

3 程序架构

根据理论模型的建立和算法的分析，设计控制棒下落时间计算软件的架构如图5所示。

核心计算程序采用C++编写，运行环境为Linux系统。

4 算例验证

为了验证理论算法的正确性，本文把编制程序的计算结果与岭澳二期核电站4号机组热态落棒试验实测数据进行对比，如表1所示。程序计算得出的控制棒下落过程加速度、速度、位移随时间的变化曲线如图6所示，该程序能够较准确地计算压水反应堆控制棒的下落时间，而且计算值较试验值偏大，说明程序计算值偏保守。

图4 程序流程图Fig.4 Program flow chart.

图5 程序架构图Fig.5 Architecture diagram.

表1 落棒时间计算结果和试验结果对比Table 1 Drop time comparison between calculation results and test results.

图6 控制棒下落过程加速度(a)、速度(b)和位移(c)的时程曲线Fig.6 Acceleration (a), velocity (b) and displacement (c) time history curve of control rod falling.

5 结语

根据对落棒时间下落过程的机理分析，利用流体力学的动量守恒方程，并且采用有限差分法，将复杂的动力学方程转化成容易求解的偏微分方程，利用时间的一阶向后差分方式进行求解计算。此外，采用面向对象的高级语言，对下落过程的各个受力模型进行抽象化处理，使得编制的程序更具有阅读性和扩展性。

通过程序计算结果与岭澳核电站二期机组热态落棒试验值进行比较，表明本文总结的落棒时间计算模型和计算方法具有可靠性，该程序可用于压水堆停堆时的安全评价。

1 孙磊, 于建华, 魏永涛，等. 控制棒组件落棒时间与历程计算[J]. 核动力工程, 2003, 24(1): 59-62 SUN Lei, YU Jianhua, WEI Yongtao, et al. Analysis of drop-time and course of control rod assembly[J]. Nuclear Power Engineering, 2003, 24(1): 59-62

2 窦一康, 姚伟达, 杨仁安, 等. 反应堆控制棒在事故工况下的落棒时间分析[F]. 中国核科技报告, 1998 DOU Yikang, YAO Weida, YANG Ren’an, et al. Computation of reactor control rod drop time under accident conditions[F]. China Nuclear Science and Technology Report, 1998

3 李庆扬, 王能超, 易大义. 数值分析[M]. 北京: 清华大学出版社, 2008 LI Qingyang, WANG Nengchao, YI Dayi. Numerical analysis[M]. Beijing: Tsinghua University Press, 2008

4 Lippman S B, Lajoie J, Moo B M, et al. C++ primer[M]. Beijing: Electronic Industry Press, 2013

5 陈卓如. 工程流体力学[M]. 北京: 高等教育出版社, 2004 CHEN Zhuoru. Engineering fluid mechanics[M]. Beijing: Higher Education Press, 2004

6 屈铁军, 陈厚群. 强地震作用下核电站控制棒下落时间分析[J]. 计算力学学报, 2001, 18(3): 364-370 QU Tiejun, CHEN Houqun. Control rod drop time analysis for the nuclear power plant under action of strong earthquake[J]. Chinese Journal of Computational Mechanics, 2001, 18(3): 364-370

CLCTL364

Calculation model of controlling rod drop time

LIU Xin CHEN Xianlong ZHANG Xiu YU Xiaolei
(China Nuclear Power Design Company Ltd., Shenzhen 518172, China)

Background:Reactor drive line is the "lifeline" of reactor control, which performs important safety functions.Purpose:To provide the basis for the drive line design assessment and nuclear power plant safety analysis.Methods:Through analyzing the stress of pressurized water reactor control rod assembly and control rod drive mechanism, established the mathematical model and the motion equation, and then solved by the finite difference method.Results:The program calculated value accords well with the test value of Ling’ao nuclear power plant phaseⅡ unit hot rod drop.Conclusion:The rod drop time of calculation program based on this method can be used to the safety evaluation when pressurized water reactor is shut down.

Controlling rod assembly, Drive line, Rod drop time

TL364

10.11889/j.0253-3219.2014.hjs.37.110604

刘新，女，1986年出生，2009年毕业于南华大学，现主要从事反应堆非堆芯热工水力计算软件研发，E-mail: liuxin_067@163.com

2014-04-04，

2014-05-08