核电厂疲劳监测系统管壁导热反演数值计算方法

刘　浪，孟阿军，凌　君，何大宇

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

核电厂疲劳监测系统管壁导热反演数值计算方法

刘浪，孟阿军，凌君，何大宇

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

核电厂疲劳监测系统中的导热反演计算是关键步骤，首先推演得到一维非稳态导热反演问题的数值计算公式，然后基于该公式，利用试验获取的管道外壁实测温度时程数据计算得到管道内壁温度，最后通过与试验实测管道内壁温度进行比较，验证了一维非稳态导热反演问题的数值计算方法的正确性。

疲劳监测系统；一维非稳态导热反演；数值计算

符号说明：

Φλ（i-1）——控制容积i-1导入i的热量，W

d U——控制容积i热力学能的增加，W

A——控制容积横截面面积，m2

λ——热导率，W/（m·℃）

Δx——控制容积长度，即空间离散步长，m

Δτ——时间离散步长，s

ρ——密度，kg/m3

c——比热容，J/（kg·℃）

ti

s——控制容积i在s时刻的温度，℃

t∞

s——环境流场在s时刻的温度，℃

a——热扩散率，m2/s，a=λ/（ρc）

To（i，j）——管道外壁的温度，℃

Ti（i，j）——管道内壁的温度，℃

0　引言

在核电厂老化与寿命管理的理论体系中，关键设备老化降质的监督监测结果是后续开展设备老化状态和寿命评估的基础。核电厂关键设备存在的主要老化机理包括疲劳、腐蚀、磨损、辐照脆化以及热老化等，其中疲劳是造成核电厂关键设备失效最主要的原因之一。目前，在世界范围内核电厂已经发生了一系列反应堆冷却剂管道与疲劳相关的事件，已影响到电厂安全、经济运行，严重时甚至引起非计划停堆。

通常，在核电厂的工程设计阶段，会按照规范要求对核一级管道进行疲劳分析，但是该分析结果不能真实反映管道的实际损耗，因为该分析过程中使用的设计瞬态进行计算，然而电厂实际运行瞬态与设计瞬态在数值上相比有一定的差距。为了获取实际运行瞬态数据，掌握管道的真实疲劳状态，基于设备的实际损伤状态，优化运行规程及在役检查大纲，合理挖掘关键管道疲劳设计的安全裕度，为电厂定期安全审查（PSR）或电厂延寿提供真实的数据支撑，部分核电国家开展了核电厂疲劳监测系统的研发工作［1］。

一般的管道计算软件如SYSPIPE，PIPESTRESS可以用设计瞬态数据开展疲劳分析得到使用系数，但是上述软件无法应用到疲劳监测系统，原因有两个:一方面，实际运行瞬态与设计瞬态相比在数据量上差别非常大，SYSPIPE，PIPESTRESS只能计算有限的瞬态个数；另一方面，设计瞬态描述的是工质温度变化，而疲劳监测系统实测到的是管道外壁温度变化，SYSPIPE，PIPESTRESS软件不具备导热反演计算的能力。综上所述，需要专门开发导热反演的计算方法，满足疲劳监测系统的需要［2-6］。

本文基于试验数据，在研究导热数值反演计算过程中，重点分析了传热数值模型、边界条件、内部节点和边界节点的温度差分方程、空间离散度对计算精度的影响。

1　分析流程

通过疲劳监测系统计算得到疲劳使用系数需要以下4个步骤，本文重点研究第2步工作。

第1步:通过疲劳监测系统在管道外壁安装的热电偶（如图1，2所示），实测获取管道外壁温度；

第2步:根据上一步实测管道外壁温度，导热数值反演计算得到管壁径向各点温度；

第3步:计算管壁径向温度梯度，包括平均温度、线性温度梯度和非线性温度梯度；

第4步:将温度梯度代入RCCM公式计算应力，最终求得疲劳使用系数。

图1　监测点示意

图2　监测点横截面图

本文重点研究第2步工作，基于试验（详见第3节）数据（包括实际测量的管道外壁和内壁温度时程数据），首先验证试验数据的准确性（详见第3.3节），即用ANSYS软件依据实测内壁温度计算得到外壁温度，并与实测外壁温度比较，判断试验数据是否准确；其次，验证一维非稳态导热反演问题的数值计算方法的正确性（详见第3.4节），即依据实测外壁温度数值反演计算得到内壁温度，并与实测内壁温度比较，判断计算精确度；然后分析壁厚的最小空间离散度（详见第3.5节），即用ANSYS软件计算不同空间离散度下壁厚的温度梯度，得到最小空间离散度，也就是说壁厚径向最少需要几个节点；随后用一维非稳态导热反演问题的数值计算方法类似第2步求的各节点温度（详见第3.6节），最后求得壁厚温度梯度（详见第3.7节），并与ANSYS软件计算的温度梯度比较，论证数值计算结果的正确性。

2　分析方法

2.1导热数值反演

2.1.1基本假设

核电厂被监测的热载荷主要有两种类型:热冲击和热分层，如图3，4所示。

图3　热冲击示意

图4　热分层示意

在热冲击的情况下管道周向温度没有差异，在热分层的情况下管道周向温度有差异，但是对于薄壁管道而言，监测点处的导热计算可以基本假设如下:

（1）管道轴向传热忽略；

（2）管道为薄壁轴对称模型，将分析对象简化为大平壁一维热传导问题。

2.1.2温度差分方程

（1）内部节点温度差分方程［7］。

如图5所示，对于无内热源的大平壁的一维非稳态导热，内部节点i所代表的控制容积（灰色）的热平衡可表述为:在s时刻，单位时间内从相邻控制容积i-1与i+1分别导入的热量Φλ（i-1）与Φλ（i+1）之和等于该控制容积热力学能的增加d U，即:

图5　一维非稳态导热内部节点温度差分方程

节点i的温度对时间的变化率采用向前差分，则热平衡方程为:

令FOΔ=AΔτ/Δx2，则上式可以整理为:

针对上式说明两点:

1）内部节点i在s时刻的温度，可以由该节点及相邻节点在s-1时刻的节点温度求出；

2）上式1-2FOΔ不能为负，称为一维非稳态导热内部节点温度方程的稳定性条件。

（2）边界节点温度差分方程［7］。

对于边界节点（如图6所示）的导热计算，同上述过程，可以得到边界节点i-1的热平衡方程为:

令FOΔ=aΔτ/Δx2，BiΔ=hΔx/λ，则上式可以整理为:

针对上式说明两点:

1）边界节点i在s+1时刻的温度，可以由该节点及相邻节点在s时刻的节点温度求出；

2）上式1-2BiΔFOΔ-2FOΔ不能为负，称为一维非稳态导热边界节点温度方程的稳定性条件。

图6　一维非稳态导热边界节点温度差分方程

（3）边界节点温度差分方程反演。

根据式（5）可以得到如下方程:

针对上式说明两点:

1）内部节点i在s时刻的温度，可以由相邻边界节点在s，s+1时刻的节点温度求出；

2）上式1-2BiΔFOΔ-2FOΔ不能为负，称为一维非稳态导热边界节点温度反演方程的稳定性条件。2.2温度梯度计算［8］

对管道进行一维非稳态导热分析后得到管壁径向的温度分布（见图7），然后根据RCC-M规范对管壁径向温度分布结果进行处理以得到温度梯度分布，即求得沿管壁径向平均温度分布T（i，j）、线性温度分布ΔT1（i，j）以及非线性温度分布ΔT2（i，j）（其中i，j分别表示载荷i与j）。本文计算温度梯度的目的是分析在一维非稳态导热数值计算过程中空间离散尺度对疲劳计算结果的影响。

图7　管壁径向温度分布

（1）平均温度T（i，j）。

（2）线性分布的温度ΔT1（i，j）。

（3）非线性分布的温度ΔT2（i，j）。

3　分析实例

3.1试验模型介绍

试验模型如图8所示。试验过程:第1步，001VP，002VP关闭，001BA内流质加热至75℃，同时003VP，004VP打开，002PO运转，用002BA中的常温流质冷却监测点管段至常温；第2步，002PO停运，003VP，004VP关闭，001VP，002VP打开，001PO运转，001BA中的加热流质进入监测点管段，形成热冲击，观察探测到的管道外壁温度数据，待温度监测数据进入稳态后，001PO停运，计算并分析监测结果。

图8　试验模型

图9　监测点横截面图

如图9所示，监测点处沿管壁径向，在管道内壁、外壁设置（焊接）热电偶，测量管壁温度变化；管道带保温，管壁厚9mm，管道材料Z2CN18.10，材料性能数据见表1。

表1　材料性能

3.2试验监测数据

如上所述，试验可以监测到管道内壁、外壁的温度数据，如图10所示。以下将根据一组试验监测到的管道外壁温度数据，进行一维非稳态导热反演的数值计算，算出管道内壁的温度，然后与试验监测的数据（如图10所示）进行对比，论证计算方法的正确性。但是在开展反演计算前，首先需要分析监测数据（如图10所示）的准确度（详见第3.3节），因为监测点处，焊接于管道外壁和内壁的热电偶会因为焊接质量问题导致监测数据失真。

图10　试验监测温度数据

3.3试验监测数据准确性分析

用ANSYS软件进行热传导计算，根据管道截面实际尺寸（外径和壁厚）确定矩形截面，依据实测内壁温度计算得到外壁温度，外壁有保温层，所以边界条件按绝热处理，ANSYS软件计算模型见图11，计算结果见图12，可以发现计算的外壁温度与实测外壁温度吻合，说明试验获取的管道内壁和外壁温度数据准确［9］。

图11　有限元模型

图12　试验监测温度数据验证

3.4一维非稳态导热反演数值计算方法正确性分析

基于试验实测外壁温度，根据第2.1.2（3）节中边界点温度差分方程反演计算管道内壁温度，与试验实测内壁温度比较（如图13所示），发现最大相差小于10℃。由于温度场数据在疲劳计算中是一个过程数据，该偏差对疲劳计算结果的影响暂时还无法判断，需要求出壁厚温度梯度后再作判断。根据第3.7节壁厚温度梯度结果可以认为，一维非稳态导热反演数值计算方法满足工程需要。

图13　反演计算内壁温度

3.5数值计算最小空间离散度分析

壁厚空间离散度对壁厚温度梯度（直接影响疲劳计算结果）有影响，本文用ANSYS软件计算不同离散尺度下壁厚温度数据，进而数值计算出温度梯度，发现平均温度、线性温度梯度受离散尺度影响较小，非线性温度梯度受离散尺度影响较大，分析可知（见图14，15），试验用的管壁至少9等分可以保证计算精度。

图14　ANSYS计算管壁径向线性温度梯度

图15　ANSYS计算管壁径向非线性温度梯度

3.6一维非稳态导热反演数值计算壁厚温度场

根据第3.5节可知，壁厚方向至少要有10个节点的温度数据，才能保证计算精度。但是根据第2.1.2（3）节中边界点温度差分方程反演的稳定性条件，时间离散尺度1 s（温度检测装置间隔1 s采集数据）的情况下，空间离散尺度最小3mm。于是管壁9等分后（如图16所示），点1，2，3，4，5，6，7的温度可以直接通过点10的温度导热数值反演计算得到，点8，9与管道外壁点10间距都小于3mm，不满足反演方程的稳定性条件，无法直接通过点10的温度导热反演计算，但是点8，9的温度既可以通过其他点反演计算得到，也可以通过线性插值简化处理，因为点7，8， 9，10温度十分接近。最终管壁上10个点温度全部得到。

图16　管壁径向离散示意

3.7壁厚温度梯度计算

利用第3.6节中一维非稳态导热反演计算得到的各节点温度，根据式（7）数值计算得到管壁径向平均温度（见图17），与ANSYS软件计算结果相比一致（ANSYS软件用试验监测的管道内壁温度数据计算得到各层温度，然后求得管壁径向平均温度）。

图17　管壁径向平均温度

利用第3.6节中一维非稳态导热反演计算得到的各节点温度，根据式（8）数值计算得到管壁径向线性温度梯度（见图18），与ANSYS软件计算结果相差小于5℃（ANSYS软件用试验监测的管道内壁温度数据计算得到各层温度，然后求得管壁径向线性温度梯度）。

利用第3.6节中一维非稳态导热反演计算得到的各节点温度，根据式（9）数值计算得到管壁径向非线性温度梯度（见图19），与ANSYS软件计算结果相差小于2℃（ANSYS软件用试验监测的管道内壁温度数据计算得到各层温度，然后求得管壁径向非线性温度梯度）。

图18　管壁径向线性温度梯度

图19　管壁径向非线性温度梯度

4　结语

（1）核电厂疲劳监测系统所涉及的管道，在热冲击和热分层两种载荷下，可以简化用一维数值模型；一维非稳态数值反演计算方法可以得到管壁各节点温度值，计算精度满足工程要求。

（2）对于一个具体的监测点，在进行导热数值反演计算前，首先需要分析最小空间离散度，目的是确保疲劳计算结果的精度，本文试验中用到的9mm厚管道，经分析至少需要10个节点。

［1］李岗，梁兵兵，马志才.稳压器波动管温度监测及寿命评估［J］.压力容器，2012，29（6）:63-67.

［2］唐永进.压力管道应力分析［M］.北京:中国石化出版社，2003.

［3］石媛媛，王明毓，汤凤，等.某研发项目稳压器排放管线应力分析［J］.压力容器，2015，32（7）:47-51.

［4］［美］Suresh S.材料的疲劳［M］.北京:国防工业出版社，1993.

［5］署恒木，黄朝琴，李翠伟.瞬态热传导问题的加权最小二乘无网格法［J］.计算力学学报，2008，25（6）:904-908.

［6］张征明，王敏稚，何树延.HTR-10核安全一级管道的力学分析［J］.清华大学学报（自然科学版），2000，40（12）:14-17.

［7］张学学，李桂馥.热工基础［M］.北京:高等教育出版社，2009.

［8］AFCEN，RCC-M，压水堆核岛机械设备设计和建造规则（2000版+2002补遗）［S］.

［9］署恒木.轴对称热传导有限元格式［J］.计算力学学报，2000，17（1）:123-126.

A Numerical Calculation of Inverse Heat Conduction for Fatigue M onitoring System of Nuclear Power Plant

LIU Lang，MENG A-jun，LING Jun，HE Da-yu
（China Nuclear Power Design Co.，Ltd.（Shenzhen），Shenzhen 518172，China）

The unsteady inverse heat conduction calculation is very important for the fatigue monitoring system，and is the focus of the thesis.The study is based on a test，the temperature of inner wall and outer wall ismeasured directly.First，the equation of numerical calculation of one dimensional unsteady inverse heat conduction is studied.Second，with the help of the equation，the temperature of innerwall is calculated by themeasured temperature of outer wall.Last，the calculated and themeasured temperature of inner wall is compared.The numerical calculation of one dimensional unsteady inverse heat conduction could be applied to the fatiguemonitoring system.

fatiguemonitoring system；one dimensional unsteady inverse heat conduction；numerical calculation

TH123；TL353

A

1001-4837（2015）11-0020-07

10.3969/j.issn.1001-4837.2015.11.004

2015-09-14

2015-10-27

刘浪（1983-），男，工程师，主要从事工程力学方面的工作，通信地址:518172广东省深圳市龙岗区黄阁北路天安数码城5号楼深圳中广核工程设计有限公司，E-mail:liulang@cgnpc.com.cn