多用途重水研究堆数字化功率调节系统仿真分析

贾玉文，毛 欢，殷浩哲，唐国静，陈 莉

(1.中国原子能科学研究院 反应堆工程研究设计所，北京 102413；2.环境保护部 核与辐射安全中心，北京 100082；3.南华大学，湖南 衡阳 421001)

20世纪90年代中国援建某国一座低温常压重水慢化和冷却、石墨作反射层、额定功率为15 MW的多用途研究堆。近期对方有意向对其进行升级改造，并明确提出要将原模拟电路的反应堆功率调节系统升级改造为全数字化控制系统。

尽管人们在将各种先进控制理论应用到反应堆功率调节系统方面进行了很多研究，如现代控制理论[1]、最优控制理论[2]、双线性控制[3]、模糊控制[4]等，但实际很少应用。截至目前，应用最为广泛的仍然是经典的PID控制。

本文基于MATLAB/Simulink平台对该多用途重水研究堆的数字化功率调节系统进行仿真研究。

1 系统组成

数字化反应堆功率调节系统采用数字计算机作为控制器，和反应堆、功率测量装置、控制棒驱动机构、调节棒等构成闭环系统。

2 系统各模型的描述

本文在MATLAB/Simulink平台下进行理论计算和仿真分析。首先建立系统内各环节的模型。

2.1 反应堆堆芯模型

1) 堆芯物理模型

采用点堆动力学方程描述堆芯中子密度。该堆为重水慢化和冷却的反应堆，γ射线可与氘核发生(γ，n)反应，在慢化剂中产生光致缓发中子，光致缓发中子对反应堆动力学的影响与缓发中子相同，故光致缓发中子也按缓发中子处理。根据半衰期的不同可将光致缓发中子分为9组[5]。故本文采用15组缓发中子点堆动力学方程：

i=1,…,15

(1)

其中：n为中子密度；ρ为反应性；Λ为中子代时间；βi为第i组缓发中子份额；λi为第i组缓发中子先驱核的衰变常量；Ci为第i组缓发中子先驱核浓度；t为时间。

2) 堆芯热工及温度反馈模型

该反应堆采用柱状燃料元件，燃料芯块为UO2，包壳为锆2合金，间隙充以氦气。共计72个燃料组件，每个燃料组件内共12根燃料元件。忽略堆芯径向和轴向的功率不均匀性，以平均功率的单根燃料元件棒作为对象构建堆芯热工模型。

芯块的传热守恒方程为：

hg2πrcil(Tfs-Tci)

(2)

Tf=(Tfs+Tfo)/2

(3)

(4)

其中：ρf为燃料密度；cf(Tf)为燃料比热容；Tf为燃料平均温度；rf为燃料芯块半径；l为单根燃料元件芯块的总长度；q(t)为燃料元件线热功率；hg为包壳与芯块之间氦气的传热系数[6]；rci为包壳内表面半径；Tfs为燃料芯块的表面温度；Tci为包壳内表面温度；Tfo为燃料芯块中心温度；df为燃料芯块直径。

包壳的传热守恒方程为：

hg2πrcil(Tfs-Tci)-h2πrcol(Tco-Tw)

(5)

Tc=(Tci+Tco)/2

(6)

(7)

其中：ρc为包壳密度；cc(Tc)为包壳比热容；Tc为包壳平均温度；rco为包壳外表面半径；h为包壳外表面与冷却剂之间的传热系数[5]；Tco为包壳外表面温度；Tw为冷却剂平均温度；kw为冷却剂导热系数；De为水力当量直径；Re为雷诺数；Pr为普朗特数。

冷却剂的传热守恒方程为：

(8)

Tw=(Twout+Twin)/2

(9)

其中：ρw为冷却剂密度；cw(Tw)为冷却剂比热容；awl为单根燃料元件棒相对应的冷却剂体积；Gw为冷却剂质量流量；Twout为冷却剂出口温度；Twin为冷却剂进口温度。

该堆温度反应性反馈主要考虑燃料的多普勒效应和慢化剂的温度效应。

ρT=ρTf+ρTw

(10)

ρTf=αf(Tf-Tf0)

(11)

ρTw=αw(Tw-Tw0)

(12)

其中：ρT为总的温度反应性反馈；ρTf为燃料温度反应性反馈；ρTw为慢化剂温度反应性反馈；αf为燃料反应性温度系数；Tf0为初始燃料平均温度；αw为慢化剂反应性温度系数；Tw0为初始慢化剂平均温度。

由于该堆反应性温度系数为负，所以即使功率调节系统不干预，反应堆也能靠温度负反馈补偿反应性扰动而达到自稳定。如当反应堆以5 MW的功率稳定运行时，分别引入1×10-4、3×10-4、5×10-4Δk/k的阶跃反应性扰动，功率响应示于图1。图2为引入5×10-4Δk/k的阶跃反应性扰动时的温度负反馈补偿扰动的曲线。

图1 不同反应性扰动的功率响应

图2 反应性扰动为5×10-4 Δk/k时的反应性变化

2.2 功率测量装置模型

功率测量装置包括电离室及相应的放大电路等。放大电路的核心为阻容电路，可等效为一阶惯性环节。

(13)

其中：n(s)为堆中子密度实际值；nm(s)为堆中子密度测量值；Tm为测量装置的时间常数；s为复变量。

2.3 控制器模型

数字化升级改造之后，采用的数字化控制系统的控制器为数字计算机，完成A/D转换及采样、定值比较、PID算法运算、D/A转换输出等功能。

1) A/D模块

A/D模块除完成A/D转换之外，还实现采样保持的功能。采样保持是将连续系统转化为离散系统的关键环节，是数字化系统与模拟系统最根本的区别。采样保持功能由零阶保持器模拟，其传递函数为：

(14)

其中，T为采样周期。

计算机运算后的控制信号经D/A模块将数字量转换为模拟输出量控制电机。

2) 定值比较器

定值比较器用于计算功率相对偏差。由于功率正比于中子密度，因此可得：

(15)

其中：ΔP为相对功率偏差；Pset为功率定值；Pm为功率测量值；nset为功率定值所对应的堆中子密度。

3) PID控制器

PID控制器接受功率偏差信号，输出棒速信号到控制棒驱动机构，其传递函数为：

(16)

其中：v(s)为棒速信号；KP为比例系数；KI为积分系数；KD为微分系数。

4) 控制器死区

为避免调节棒频繁动作，控制器输出设置有死区。在相对功率偏差为±0.1%的范围内，控制器输出的棒速信号为零。

2.4 控制棒驱动机构模型

控制棒驱动机构执行部件为步进电机，通过齿轮减速箱减速，带动绳轮旋转，由钢丝绳将旋转运动变为调节棒的上下移动。从速度到位移可等效为积分环节。

(17)

其中，L(s)为调节棒位移。

根据步进电机的步距角和减速箱的减速比，计算得到步进电机每转动一步对应的调节棒位移为0.025 mm。调节棒的上提最大速度为60 mm/s，下插最大速度为500 mm/s。

2.5 调节棒反应性模型

堆内有两根调节棒，一根参与调节反应性，另一根备用。根据物理数据采用6次多项式对单根调节棒在全行程(1 360 mm)内的反应性价值进行拟合：

(18)

其中：l为调节棒的位置；ai为多项式系数。

根据拟合公式得到的调节棒价值曲线如图3所示。

图3 调节棒价值曲线

3 系统稳定范围理论计算

将系统进行适当简化：1) 暂不考虑棒速饱和、控制器死区等非线性因素；2) 调节棒模型采用线性段微分价值代替；3) 暂不考虑温度反馈。对简化的连续系统模型进行离散化，采用离散控制系统稳定性判据“保持离散系统稳定的充分必要条件是闭环离散系统特征方程的根的模均小于1”[7]，计算出不同采样周期T下的比例系数，结果列于表1。

表1 不同采样周期下的比例系数最大值

最佳控制器参数的整定只需在相应范围内进行选择，可有效避免参数整定的盲目性。

4 仿真分析

4.1 关于PID控制器的进一步说明

本方案采用的PID控制器仅采用比例控制，原因如下。

从功率偏差到反应性扰动输入点之间存在一个积分环节(从速度到调节棒的位移)，使得系统对阶跃反应性扰动的稳态误差为零[7]，因此没有必要增加积分控制。而且增加积分环节会增加系统阶数，可能会带来系统不稳定的因素。

微分控制很容易放大高频噪声，降低系统的信噪比，使系统抗干扰能力下降。工程实际应用中，也要求控制系统尽量简单可靠。

采样时间和PID控制器参数的整定通过Simulink仿真分析完成。离散控制系统的Simulink模型如图4所示。

4.2 稳定性分析

假设反应堆功率为满功率(15 MW)，在1 s时引入5×10-4Δk/k的阶跃反应性扰动，KP=2 000时的功率响应曲线如图5所示。冷却剂出口温度稳定性分析如图6所示。

由图5、6可见，功率超调量约为6%定值功率，稳定时间小于2 s，冷却剂出口温度也随时间逐渐稳定下来。因此闭环离散控制系统是稳定的。

4.3 采样周期和控制器参数的整定

1) 不同采样周期下最佳KP的整定

假设反应堆功率为满功率(15 MW)，在1 s时引入5×10-4Δk/k的阶跃反应性扰动，不同采样周期下采用不同比例系数的功率响应如图7所示。

图4 数字化功率调节系统Simulink模型

图5 功率响应曲线

图6 冷却剂出口温度稳定性分析

采样周期：a——10 ms；b——100 ms；c——500 ms；d——1 s

由图7可知：当KP较大时，虽然超调量较小，但振荡较为厉害，且稳定时间较长；当KP较小时，虽然振荡可得到有效控制，但超调量和稳定时间皆变差。因此综合考虑超调量、振荡和稳定时间3项性能指标，保持系统性能最佳的参数可做如下选择：T=10 ms时，KP=3 500；T=100 ms时，KP=2 000；T=500 ms时，KP=1 000；T=1 s时，KP=500。

2) 不同采样周期下控制器性能的比较

对采用不同采样周期(10、100、500 ms)的控制器的性能进行比较，结果如图8所示。由图8可见：每个控制器的控制器参数均为其相应的最佳参数；减小采样周期有利于改善系统的控制性能，但采样周期过短，对系统控制性能的提高不是很大，反而会显著增加控制系统的计算负荷，采样周期过长，控制性能变差，甚至无法满足要求。

综合上述分析，反应堆数字化功率调节系统PID控制器的采样周期定为100 ms，KP定为2 000。

4.4 斜坡反应性扰动调节特性

假设反应堆功率为满功率(15 MW)，在10 s时引入速率为2×10-4(Δk/k)/s的反应性，持续10 s，反应堆功率、燃料平均温度、冷却剂平均温度的变化如图9所示。

图8 不同采样周期的控制器性能比较

由图9可见，在该功率调节系统的控制作用下，引入斜坡反应性时，功率始终维持在一定的(101%满功率)范围内。当斜坡反应性作用消失后，功率很快就被控制在功率定值附近。功率稳定后，燃料平均温度和冷却剂平均温度也先后趋于稳定。

4.5 随动调节特性

假设反应堆初始功率为10 MW，功率定值以10 kW/s的速率增加到11 MW，反应堆功率、燃料平均温度、冷却剂平均温度的变化如图10所示。

图9 斜坡反应性扰动下功率、燃料平均温度和冷却剂平均温度的响应

图10 随动调节时功率、燃料平均温度和冷却剂平均温度的响应

由图10可见，在该功率调节系统的控制作用下，反应堆功率可快速平稳地跟随功率定值变化，并可将反应堆功率维持在最终的功率定值附近。

5 结论

1) 采用PID控制器可满足该反应堆数字化功率调节系统控制性能的要求；

2) 采用适当的采样周期和PID控制器参数可克服反应性扰动的影响，保持系统稳定运行，并可快速平稳地跟随功率定值变化。

理论计算的结果对于控制参数的整定具有一定的指导意义。

参考文献：

[1] 杨志军，刘隆祉，胡桂芬. 核供热堆功率自动控制系统仿真研究[J]. 原子能科学技术，2000，34(增刊)：94-97.

YANG Zhijun, LIU Longzhi, HU Guifen. The simulation study on the nuclear heating reactor’s power auto-control system[J]. Atomic Energy Science and Technology, 2000, 34 (Suppl.): 94-97(in Chinese).

[2] 赵福宇，曹艳. 反应堆最优功率控制系统的设计[J]. 核科学与工程，1999，19(1)：39-43.

ZHAO Fuyu, CAO Yan. Design of reactor power optimized control system[J]. Chinese Journal of Nuclear Science and Engineering, 1999, 19(1): 39-43(in Chinese).

[3] 邓琛，张琴舜. 双线性控制技术在反应堆功率调节系统中的应用[J]. 上海交通大学学报，1997，31(9)：149-151.

DENG Chen, ZHANG Qingshun. Application of dual-linear control technology in reactor power control system[J]. Journal of Shanghai Jiao Tong University, 1997, 31(9): 149-151(in Chinese).

[4] RUAN Da, LI Xiaozhong. Comparative study of fuzzy control, PID control, and advanced fuzzy control for simulating a nuclear reactor operation[J]. International Journal of General Systems, 2000, 29: 263-279.

[5] 仲言. 重水研究堆[M]. 北京：原子能出版社，1989.

[6] LANNING D D, HANN C R, GILBERT E S. Statistical analysis and modeling of gap conductance data for reactor fuel rods containing UO2pellets, BNWL-1832[R]. USA: Battelle Pacific Northwest Laboratories, 1974.

[7] 徐薇莉，曹柱中，田作华. 自动控制理论与设计[M]. 上海：上海交通大学出版社，2001.