基于薄板样条函数插值的核电站堆芯温度场软测量方法与流程

本发明属于热工测量技术领域，特别涉及基于薄板样条函数插值的核电站堆芯温度场软测量方法。

背景技术：

世界第一座试验性核电站于1954年6月在前苏联奥布宁斯克建成，它的建成标志着核电时代的到来。在随后的二十年，核电发展逐渐步入高潮，许多先进的工业国家建了一批核电站。核能逐渐成为我国及全世界能源利用的一个重要组成部分。

但是随着核能的发展，核电事故的接连发生却给人们带来了巨大的恐慌，核安全问题受到人们的关注。核反应堆堆芯是核电站的要害部位，当堆芯温度超过限值时会发生熔融，造成核泄漏。在任何情况下保证堆芯产生的热量能及时地输出是保证核电站安全的关键。只有堆芯内的燃料元件处于适当的冷却状态，热量被及时带出，才能确保反应堆及核电站的安全。

国内外对于温度场重构方法的研究，大多数解决的是炉膛流动介质的温度场重构问题，但炉膛结构和反应堆堆芯结构存在较大差异，并且温度的测量方法不同等因素，这些方法有借鉴之处，但无法直接应用。另外，反应堆堆芯温度测点的位置和数量受到实际条件的限制，只能布置在冷却剂出口处，无法监测堆芯组件的温度，不能全局观测堆芯温度分布，存在监视盲区。

薄板样条函数是一种通过空间上所有散点找到一个最小弯曲的光滑曲面的插值方法，可以通过三个空间散点得到平面。因此，如何通过薄板样条函数插值的方法，利用现有的传感器测点数据建立反应堆堆芯温度场软测量模型，实现对堆芯内部任意位置温度值的估计，对保障核电站运行的安全性具有实际意义的。

技术实现要素：

发明目的：针对现有技术中存在的问题，本发明提供基于薄板样条函数插值的核电站堆芯温度场软测量方法，该方法能够以布置在堆芯冷却剂出口的传感器采集到的温度为初始数据，在现有压水堆堆芯功率分布模型的基础上，结合稳态传热热量守恒方程，应用薄板样条函数(Thin Plate Spline)插值法，建立堆芯温度场重构的物理模型和数学模型，最终实现对压水堆堆芯冷却剂和堆芯燃料组件的三维温度分布的软测量结果。

技术方案：为解决上述技术问题，本发明提供基于薄板样条函数插值的核电站堆芯温度场软测量方法，包括以下步骤：

步骤一：如图2所示，基于核电站反应堆结构及堆芯冷却剂出口处的温度传感器测点的分布特点建立坐标系(参见图2)，将传感器测点温度作为离散点，采用薄板样条函数(Thin Plate Splines，TPS)插值的方法对测点温度进行数值计算处理，对其他未布置温度传感器的冷却剂出口温度进行估计，初步重构得到堆芯冷却剂出口整个截面的温度场的分布；

步骤二：通过反应堆结构及其回路系统中冷却剂的流动方式(参见图3)、热量传输特性和燃料元件的传导特性，建立堆芯温度计算模型和堆芯热传导模型，并计算得到堆芯的径向功率Q_ij、堆芯的轴向功率Q_k和每个节块的功率Q_ijk；

步骤三：根据各通道的冷却剂进/出口温度，取其平均值调用冷却剂热力参数(比热容、密度)，采用换热公式Q_ij＝c_ij·w_ij·ρ_ij·ΔT_ij，计算出冷却剂在各通道的流量w_ij；式中,ΔT_ij是各通道冷却剂进、出口温度差，Q_ij是堆芯径向功率，c_ij是各通道冷却剂的比热容，ρ_ij是各通道冷却剂的密度；如图2中，冷却剂的各通道流量即为箭头所指向的冷却剂在压水堆中流动的所有通道流量。

步骤四：利用步骤二得到的每个节块的功率Q_ijk和步骤三所得到的各个通道的冷却剂流量值w_ij，从入口或出口开始，将上一个节块的温度及对应的冷却剂密度和比热带入换热公式Q_ijk＝c_ijk·w_ij·ρ_ijk·ΔT_ijk中中计算下一个节块温度的初值，再根据温度初值修正本节块冷却剂的密度和比热容，从而得到进一步精确的节块温度，如此反复迭代逼近直至下个节块的温度满足精度要求。以此类推，最终得到整个堆芯的冷却剂温度分布；式中，c_ijk为每个节块冷却剂的比热容，w_ij为每个冷却剂通道的流量，ρ_ijk为每个节块冷却剂的密度，ΔT_ijk为每个节块冷却剂的进出口温差。

步骤五：将步骤四所得到的核电站堆芯冷却剂的三维温度带入步骤二中的堆芯传导模型中，计算得到堆芯燃料组件各部分的三维温度分布。

进一步的，步骤一中采用薄板样条函数插值的方法对测点温度进行数值计算处理的具体方法如下：

基于核电站反应堆结构及堆芯冷却剂出口处的温度传感器的分布特点建立坐标系(参见图2)。其中传感器分布坐标的原点为堆芯流量所有通道的中心通道；具体参照图2所示，堆芯里对称布置着很多燃料组件，这些组建彼此之间都有通道，所有通道的最中心极为原点，根据原点建立x轴和y轴。核电站冷却剂出口测量得到的n个温度场离散点数据为x_i,y_i,z_i，其中i＝1，2…n，x、y为传感器分布坐标，z为对应坐标的温度。薄板样条函数插值得到的所需要重构的堆芯冷却剂出口截面的二维温度分布函数表示为：式中，w_i为权值，U为TPS核函数，其表达式表示为：U(d)＝d²log(d)；式中，d_i＝(x-x_i)²+(y-y_i)²。函数表达式中的各个系数w_i、a₁、a₂、a₃由方程组z(x_i,y_i)＝z_i(i＝1…30)，计算得出；之后可以根据任意坐标位置的x，y值按照函数表达式计算出对应的出口温度z的值。

进一步的，步骤二中的堆芯温度计算模型由以下几个部分组成：

A.堆芯通道模型：堆芯内冷却剂的横向流动较纵向流动小得多，可以不考虑通道间的质量、动量和能量交换，将堆芯内冷却剂的流动通道看作是单通道模型；

B.堆芯节块划分(参见图4)：使模型变为有限元；径向划分时，1个节块对应1个燃料组件；轴向平均划分时，尽可能使节块形状接近立方体以保证在各方向采用相同的计算方法，划分为16个节块；

C.功率分布模型：根据全堆芯各位置的中子功率预测值得到堆芯归一化的径向功率和轴向功率；每一节块的功率通过下式计算得到：

式中Q_ijk为每个节块的功率，Q为整个堆芯的已知功率，Q_ij为堆芯的径向功率，Q_k为堆芯的轴向功率，g为堆芯的热效率；

进一步的，所述步骤二中的堆芯热传导模型(参见图5)由以下几部分组成：

A.堆芯燃料元件传热模型：圆柱坐标系下轴向无限大一维热传导模型求解得到4πλ_rl(T₁-T₂)＝Q；

B.燃料元件与包壳之间间隙的传热模型：(T₂-T₃)h_gA₁＝Q；

C.包壳管内的传热模型和包壳外表面：圆柱坐标系下无内热源一维热传导模型求解得到

D.冷却剂的对流换热模型：视为均匀壁温边界条件的大空间自然对流(T₄-T_avg)h_fA₂＝Q。

式中，T₁为燃料芯块的中心温度；T₂为燃料芯块的边缘温度；λ_r为燃料芯块的热导率；l为燃料棒的高度；T₃为包壳内表面温度；h_g为燃料芯块与包壳之间的等效传热系数；A₁为芯块外表面传热面积；T₄为包壳外表面温度；λ_c为包壳材料的热导率；r_s为包壳的外径；r_i为包壳的内径；T_avg为一回路冷却剂平均温度；h_f为包壳外表面与冷却剂之间的对流换热系数；A₂为包壳外表面传热面积。

与现有技术相比，本发明的优点在于：

本发明基于计算机信息处理技术进行堆芯温度场重构，全景化地展示堆芯的温度分布，为核反应堆安全运行提供可视化的温度场图像，同时为开发核反应堆燃料元件(或燃料棒)的寿命损耗管理提供科学依据，更进一步地保障了堆芯运行的安全性和经济性。因此，本发明在未来核电站的堆芯温度检测和安全运行方面必将具有广泛的应用前景，产生显著的社会效益。且本发明还具有如下的特点及优点：(1)克服了在温度测点有限、已知数据较少的情况下，无法直接选择合适的离散数据处理方法的困难；(2)相对于其他插值/拟合方法，该方法得到冷却剂出口温度的测量值与计算值插值的平均值E和标准差σ最小；(3)本发明且能够利用堆芯出口的点温度测量值对计算模型进行自适应修正，减小了数据误差，很好地表征了堆芯内部的温度分布状况。

附图说明

图1为本发明的总体流程图；

图2为压水堆堆芯温度测点布置图；

图3为压水堆冷却剂流动路径；

图4为堆芯的三维节块划分图；

图5为堆芯热传导示意图；

图6为薄板样条插值方法计算流程图。

具体实施方式

下面结合附图和具体实施方式，进一步阐明本发明。

本专利实现方案的主要流程具体如下(参见图1)：

(1)采用的薄板样条函数插值的方法对n个测点温度进行计算处理的具体的实现如下(参见图6)：

核电站冷却剂出口测量得到的n个温度场离散点数据为x_i,y_i,z_i，其中i＝1，2…n，x、y为传感器分布坐标，z为对应坐标的温度。堆芯冷却剂出口截面的二维温度分布函数表示为：矩阵C中存储n个空间散点的温度数据：

引入正则化参数λ，以解决未知的薄板样条函数权值w以及线性方程组。如果λ为0就是普通的薄板样条函数(TPS)，如果λ越来越大，TPS就退化为均方误差最小平面(平面的弯曲量为0)。

式中的K为存储TPS核函数的矩阵，P和O是辅助子矩阵，w为存储待求参数权值w_i的列向量，a为存储待求参数a₁、a₂、a₃的列向量，v为存储堆芯出口测点温度的列向量，o是辅助列向量，具体求解如下：

K_ij＝U{|(c_i1c_i2)-(c_j1c_j2)|}+I_ij·α²·λi,j∈[1…n],λ>0

|(c_i1c_i2)-(c_j1c_j2)|表示欧几里德范数，表示(c_i1，c_i2)和(c_j1，c_j2)的直线距离；

当知道w和a的值后，可以根据任意的x，y值按照公式计算出z的值，即可得到任意出口位置的冷却剂温度值。

(2)调用冷却剂热力参数，对每个冷却剂通道采用换热公式Q＝c·w·ρ·ΔT得到相应的冷却剂流量。由于涉及调用冷却剂的比热容c、冷却剂的密度ρ等冷却剂的热力参数，采用由Prof.Hans-Joachim Kretzschmar编写的用于调取水和蒸汽物性参数的matlab程序FluidLAB LiblF97，以保证在计算过程中使用该程序进行冷却剂热力性质计算时能达到计算精度。具体调用方式为：

冷却剂比热容调用函数：cpwaterCE1＝cp_pTx_97(p,T,X)

冷却剂密度调用函数：rwaterCE1＝1/v_pTx_97(p,T,X)

式中，p为冷却剂压力，T为冷却剂温度，X为冷却剂质量。

(3)再得到了每个通道的冷却剂流量后，从冷却剂入口温度纵向第1层节块开始正向逐层计算相应节块的出口温度，或者从重构得到的冷却剂出口温度纵向第16层的节块开始反向逐层计算相应节块的入口温度。主要的数值计算方法是迭代法，其基本思想是逐次逼近，先取一个粗糙的近似值，然后用同一个递推公式，反复校正此初值，直至达到预定精度要求为止。

本步骤中，需要根据上一个节块的温度及对应的冷却剂密度、比热容计算下一个节块的温度的初值，再根据计算的温度初值修正本节块冷却剂的密度和比热容，然后得到进一步精确的节块温度，如此反复迭代逼近直至下个节块的温度满足精度要求。，以此类推最终得到整个堆芯的冷却剂温度分布。此时，所采用的换热公式相应地变为Q_ijk＝c_ijk·w_ij·ρ_ijk·ΔT_ijk，式中Q_ijk为每个节块堆芯燃料组件的功率，c_ijk为每个节块冷却剂的比热容，w_ij为每个冷却剂通道的流量，ρ_ijk为每个节块冷却剂的密度，ΔT_ijk为每个节块冷却剂的进出口温差。

(4)对堆芯传导模型求解得到：4πλ_rl(T₁-T₂)＝Q，(T₂-T₃)h_gA₁＝Q，(T₄-T_avg)h_fA₂＝Q。式中，T₁为燃料芯块的中心温度；T₂为燃料芯块的边缘温度；l为燃料棒的高度；λ_c为包壳材料的热导率；r_s为包壳的外径；r_i为包壳的内径。将上步所得的每个节块冷却剂的温度带入T_avg中，就可进一步计算得到堆芯燃料组件各部分的温度。

以上所述仅为本发明的实施例子而已，并不用于限制本发明。凡在本发明的原则之内，所作的等同替换，均应包含在本发明的保护范围之内。本发明未作详细阐述的内容属于本专业领域技术人员公知的已有技术。