基于子通道模型的超临界二氧化碳燃料组件热工特性分析

郭 嘉,赵富龙,*,宁可为,方华伟,易经纬,田瑞峰,谭思超

(1.哈尔滨工程大学 核科学与技术学院,黑龙江 哈尔滨 150001;2.中国核动力研究设计院,四川 成都 610213)

超临界二氧化碳(supercritical carbon dioxide,简称S-CO2)反应堆是近些年备受关注的四代堆型,在小型堆的研究及开发中被提及,其主要应用于布雷顿循环系统,是一种高效清洁的能源和动力转换装置,并且可以与多种形式的热源结合,实现联合循环发电[1]。相对于蒸汽朗肯循环,S-CO2布雷顿循环系统具有循环效率高、回路简单、结构紧凑等特点,在小型堆供电及动力装置方面具有良好的发展前景[2]。S-CO2兼具液体和气体的特性,非常适用于作为反应堆一回路的冷却剂,其优点包括热导率高、黏度小、无相变、以及良好的化学稳定性、易于获取且成本低廉[3-4]。

目前对S-CO2燃料组件应用子通道分析方法来进行热工水力分析的研究十分有限,研究集中在使用S-CO2作为工质的布雷顿循环系统的特性分析[5]。为了分析S-CO2工质在反应堆内的热工水力特性,获取二维尺度上的参数,本文采用适用性和可靠性高的子通道分析方法进行研究。首先选定一种燃料组件,并对其进行子通道划分,通过连续性方程、能量守恒、动量守恒方程等构建子通道模型,并使用Modelica语言开发燃料组件子通道分析程序[6-7]。子通道分析方法的收敛性、重复性和可扩展性要优于CFD计算软件,并且耗费的计算资源较少,方便后续的程序开发和程序扩展[8],为后续容积控制系统和布雷顿循环系统的程序开发和理论研究提供了技术支持[9]。子通道分析程序由于其易收敛且准确性高的特点,能够为流动换热特性的研究及堆芯瞬态运行特性分析提供重要的技术支持[10]。

目前关于S-CO2反应堆的研究集中于流动换热研究和三维CFD模拟。Sharabi等[11]针对S-CO2反应堆模型开展流动换热特性分析,通过三维CFD计算对传热恶化机理和湍流交混下冷却剂流动进行了研究。Sun等[12]研发了子通道分析程序SUBAC,主要应用于钠冷快堆带绕丝的燃料组件的分析。Ding等[13]使用Modelica语言构造了S-CO2反应堆瞬态模拟程序,进行了瞬态工况的研究并应于与布雷顿循环系统的开发。Ming等[14]开发了S-CO2为工作介质的布雷顿循环系统分析程序,研究了系统运行特性和设备之间的耦合特性。陈立新等[15]开发的子通道程序PRTHA应用在西安脉冲堆上,与实验结果比较证明了开发程序的合理性和可行性。苏开放等[16]对氦氙子通道湍流交混的研究,通过CFD计算拟合出适用于反应堆的交混系数关系式,量化分析了小型堆的交混效应。吴宗芸等[17]开发了钠冷快堆子通道分析程序SPLICA,并通过与棒束实验的对比,验证了子通道分析程序的准确性。综上所述,关于S-CO2的程序开发的研究集中于布雷顿循环系统的建立及参数提取,缺乏S-CO2反应堆的精细化模型的相关研究,为弥补相关研究的不足,本文将使用子通道分析方法对S-CO2反应堆进行程序开发。因此,本研究基于Modelica语言开发S-CO2布雷顿反应堆系统子通道分析程序BRESA,对S-CO2工质反应堆的轴向和径向的物理热工参数进行计算和分析,研究各典型通道的焓升及温升变化并验证程序的正确性。将开发的程序进行流量变化瞬态工况研究,分析系统的瞬态响应特性,得到瞬态工况下各典型通道的流量焓以及压降变化,实现子通道模型的功能,旨在为S-CO2反应堆系统瞬态运行特性研究及安全特性分析提供参考。

1 数值计算方法

1.1 子通道模型

子通道分析方法是分析燃料组件热工水力特性中有效快捷的方法之一。本文通过将燃料组件在高度上分割成多段,再根据划分的子通道将整个燃料组件中冷却剂通道划分为控制体,径向通过子通道划分,轴向上根据高度划分成10段。燃料组件根据子通道方法控制体划分如图1所示。

图1 子通道控制体划分

构建完燃料组件控制体的几何模型后,需要对每个控制体建立物理热工模型,其方法为在每个控制体上建立子通道方程,控制体对应的燃料棒附近建立反应堆释热的方程,给定初始值和边界条件,最终得到每个控制体上的物理热工参数,对控制体建立的子通道方程包括质量守恒方程、能量守恒方程、轴向和径向动量守恒方程,这些方程不仅有控制体本身的方程,还涉及到不同控制体之间的方程,最后将所有控制体的方程联立并求解,可得到整个燃料组件的物理热工参数。

由连续性方程得到每个控制体的质量守恒方程、能量守恒方程、轴向动量守恒方程、径向动量守恒方程[18-20]如下:

mi,n=mi,n-1-mi,n+1-ωin,jn+ωkn,in

(1)

qin-∑ϖin,jn(Hin,out-Hin,in)+

∑ωin,jn(Hin,in-H*)-∑Css(Tin,out-Tin,in)

(2)

(3)

(4)

式中:mi,n为单位长度的轴向流量,kg/(m·s);ωin,jn为单位长度的横向流量,kg/(m·s);uin,in为控制体的入口速度,m/s;Hin,in为控制体的入口的焓,kJ/kg;qin为线功率,kW/m;ϖin,jn为湍流交混流量,kg/(m·s);Tin,in为入口温度,K;Cs为相邻通道热导率,W/(m2·K);s为相邻通道间距,m;pin,in为控制体的入口压力,Pa;pin,out为控制体的出口压力,Pa;f为摩擦阻力系数;ft为湍流交混系数;z为控制体高度,m;d为等效直径,m;ρ为流体密度,kg/m3;uin,jn为控制体间横向流速,m/s;KG为横向形状阻力系数;s为控制体的宽度,m;下标in,in代表控制体i,n的入口,下标in,out代表控制体i,n的出口,角标*代表相邻控制体间交换的量,下标i,n代表控制体i,n,下标i,n-1代表控制体i,n-1,下标in,jn代表在控制体i,n和j,n之间进行流动。

1.2 几何模型

CO2的临界温度和临界压力分别为304.13 K和7.38 MPa,当流体温度和压力分别高于临界温度和临界压力时,被视为超临界状态,其物理性质在此状态下具有良好的稳定性,本文采用超临界态CO2作为工质,由于设计的燃料组件的温度和压力远大于临界点,所以不需考虑跨临界问题。本文设计的燃料组件为六角形燃料组建盒,燃料棒的排列方式为三角形排列,采用37棒的燃料组件,根据燃料组件对称性,对1/6堆芯进行作图,并对冷却剂通道进行逐一编号,包括10根燃料棒,并划分成14个子通道,如图2所示。

图2 1/6燃料组件子通道编号图

其中燃料组件的高度为1.0 m,燃料组件的对角线距离为164.4 mm,对边距离为142.4 mm,燃料棒直径为20 mm,燃料棒之间的间距为2.6 mm,燃料组件的P/D=1.13,子通道1～9的截面积为64.09 mm2,子通道11、12、13的截面积为127.68 mm2,子通道10、14的截面积为37.43 mm2。即子通道流通面积11、12、13最大,1～9次之,10、14流通面积最小。几何参数列于表1。

表1 子通道的几何参数

1.3 计算方法

本文基于Modelica语言开发了S-CO2布雷顿反应堆系统子通道分析程序BRESA,首先根据37棒模型,取1/6堆芯,对子通道进行编号,沿轴向将每个子通道划分为10个控制体,合理分配划分控制体的数量,在节省计算资源的情况下保证计算精度,计算流程如图3所示。

图3 计算流程

计算步骤如下:1) 首先设定初始值,确定边界条件;2) 需要先进行物理计算,求出燃料组件的释热率、冷却剂的吸热量等参数;3) 根据冷却剂吸热量求解子通道控制体的能量、动量、质量守恒方程,联立求解;4) 对所有控制体的方程进行求解,求出第1个时间步长整个燃料组件的参数;5) 每个时间步长进行1次迭代计算,将上一时间步长的值代入到下一时间步长中;6) 迭代计算直至稳定,即燃料组件中物理计算的释热量和热工计算中冷却剂的吸热量在空间上达到稳定,输出结果,实现子通道程序功能。

2 稳态特性分析

2.1 子通道出口参数分析

对子通道分析程序BRESA设定入口边界条件如下:冷却剂入口焓为795.73 kJ/kg,燃料组件的入口流量为5.553 kg/s,入口温度为628.15 K,入口压力为19.9 MPa,入口处冷却剂密度170.86 kg/m3,燃料组件功率为1 000 kW,出口边界条件为焓976.00 kJ/kg、温度775.15 K、压力19.85 MPa。

稳态计算结果的可靠性和准确性已经在前期的论文中经过验证[20],受篇幅限制不做过多介绍。

在计算达到稳定后,子通道的流量也趋于稳定,冷却剂各通道进出口流量如图4所示。其中,相同截面积的子通道,例如通道1,2,3,4…,9的进口流量相同,各子通道流通面积如表1所列,子通道流通面积与本节中子通道流量呈正相关,所以入口流量只与子通道的流通面积相关。冷却剂经过反应堆内加热和交混之后,子通道出口流量相对入口流量发生了波动,子通道的出口流量相对进口流量的变化量约为0.2%～1%,证明了各子通道之间的流量交混效应和横向流动现象的存在,且会对燃料组件不同通道之间的能量交换造成一定的影响。

图4 冷却剂子通道入口、出口流量

冷却剂各子通道出口焓如图5所示,各通道入口焓均为795.42 kJ/kg。观察可知,越靠近燃料组件中心处的子通道出口焓越高,即冷却剂在经过燃料组件的加热后的焓升越大,其原因为靠近堆芯中心处的子通道相邻的燃料棒多,因此从燃料棒吸收的热量较多,而堆芯边缘处子通道相邻的燃料棒少,吸收的热量较少,这一分布特性符合燃料组件出口处温度和焓的分布规律。

图5 冷却剂子通道出口焓

冷却剂各子通道入口、出口密度如图6所示。各通道入口密度均为170.86 kg/m3,出口密度分布趋势为越靠近燃料组件中心处的子通道出口密度越小,其中通道1出口密度为124.99 kg/m3,通道14的出口密度为142.82 kg/m3。由此可知,在经过燃料组件的加热后,随着冷却剂焓的升高,冷却剂密度逐渐降低,焓高的子通道密度较小,焓低的子通道密度较大。

图6 冷却剂子通道入口、出口密度

2.2 典型通道沿轴向参数分析

子通道分析程序BRESA的典型通道沿轴向焓升变化如图7所示。由于相临近的子通道焓近似,为便于观察,取典型通道1、5、11、14进行绘图。其中子通道1的焓升最大,为233.93 kJ/kg,通道5焓升为209.44 kJ/kg,子通道11焓升为161.70 kJ/kg,通道14焓升最小,为121.60 kJ/kg。观察可知,冷却剂在轴向上焓变化呈现在进出口焓升较小、轴向中间处焓升较大的特点,即燃料组件的热量集中在轴向的中间处。各通道的吸热量也存在差异,通道1吸热最多,轴向上温度变化剧烈,靠近堆芯边缘的通道11和14吸热量小,轴向上焓变化平缓。

图7 典型通道沿轴向焓升变化

典型通道沿轴向温升变化如图8所示。其中子通道1温升最大,为189.69 K,子通道14温升最小为99.06 K。观察可知,冷却剂在轴向上的的温度分布和焓分布近似,靠近中心的冷却剂通道温升大,靠近边缘的冷却剂通道温升少。且在轴向中间处温升变化较大。

图8 典型通道沿轴向温升变化

对比图9中Zhang等的S-CO2子通道研究中各典型通道的轴向温度分布结果[21],由于S-CO2工质在该设计工况下,温度和焓的变化趋势基本一致,则可用Zhang等的温度分布结果来验证本文中温度和焓分布是否合理。本文中典型通道的冷却剂焓和温度在轴向上的变化趋势与Zhang等的结果一致,均为斜率先增加后减少的上升曲线,并且程序得到的结果符合理论分析的结果,满足中间通道温升和焓升较高,边缘通道温升和焓升较低的分布特点,证明了子通道分析程序结果的合理性,为后文中瞬态计算提供了理论基础和程序支持。

图9 冷却剂轴向温度分布[21]

3 瞬态特性分析

3.1 瞬态工况下流量变化分析

本文针对子通道燃料组件模型,确定模型的边界条件,引入流量变化的瞬态工况,观察系统参数随时间变化情况,分析系统的瞬态响应特性。

对子通道分析程序BRESA设置流量变化如下:设定额定流量为5.553 kg/s;从1 000 s到1 500 s时间段流量线性增加,从5.553 kg/s增加到11.106 kg/s;从1 500 s到2 000 s时间段流量线性减少,从11.106 kg/s减少到8.330 kg/s;之后在极短的时间内急剧减少到额定流量5.553 kg/s。共计3个瞬态过程:一次流量线性增加,一次流量线性降低,一次流量突降。观察子通道燃料组件模型参数随时间的变化。

图10为在设置的瞬态工况下,典型通道的流量变化。由于截面积相同的子通道流量近似,为便于观察,取典型通道1、5、11、14进行绘图。其中,子通道11流量最大,子通道1和通道5流量近似,子通道14流量最少。观察可知,每个子通道的流量变化趋势基本与燃料组件的流量变化趋势保持一致,在一定范围内的总流量改变,对冷却剂通道的流量分配基本不会产生的影响,且各冷却剂通道的流量也呈比例变化。

图10 典型通道流量变化

3.2 瞬态工况下焓变化分析

在设置的流量变化瞬态工况下,反应堆会进行快速响应,功率会进行快速的变化,从而使冷却剂焓发生变化,由于燃料和冷却剂的温度负反馈系数的影响以及采用的PID算法调节的功率控制系统,反应堆的功率和冷却剂焓不是随流量线性变化的。图11为在流量瞬态变化下的功率变化曲线,式(5)为PID控制方程[14]:

图11 燃料组件功率及流量变化曲线

(5)

式中:u(t)为反馈信号;e(t)为相对误差;Kp为比例增益;Tt为积分时间常数;Td为微分时间常数。

在设置瞬态工况下,冷却剂流量线性增加过程中,冷却剂温度和焓会下降,通过温度负反馈的影响,反应堆会通过功率控制系统引入正反应性,从而使反应堆的功率上升,冷却剂温度和焓又随之上升,温度负反馈又会导致功率有小幅下降,由于冷却剂流量还处于增加的过程,冷却剂温度和焓上升一段时间后又会继续降低,因此,功率呈现波动上升的趋势,冷却剂焓呈现波动下降的趋势。冷却剂流量线性减少过程同理,功率波动下降。冷却剂流量突降时,反应堆功率也出现突降,然后产生大幅波动,最终稳定到一个新的水平。

本文的控制目标为在流量变化的瞬态工况下,保证反应堆功率和冷却剂焓能够随流量变化,避免发生反应堆超功率或冷却剂温度过高等工况,确保反应堆的安全运行。

图12为子通道分析程序BRESA的典型通道的焓变化。由于相临近的子通道焓近似,为便于观察,对典型通道1、5、11、14进行绘图。观察可知,在流量变化的瞬态工况下,冷却剂通道的焓会发生一定量的波动。流量线性增加的过程中,冷却剂通道的焓呈现波动式下降趋势,这是因为发生流量降低的瞬态过程中,冷却剂流量降低,为了防止冷却剂吸热量降低导致堆芯温度上升,反应堆功率会缓慢下降,同时冷却剂通道吸热量降低,冷却剂的焓也会缓慢下降。流量线性减少的过程中,同理,冷却剂通道的焓会波动式上升。

图12 典型通道焓变化

但当冷却剂流量急剧下降时,冷却剂焓会急剧升高,在流量下降瞬间的典型通道焓升达到90～120 kJ/kg,焓波动幅最大达到12%,其原因为冷却剂流量的急剧减少,会使冷却剂的吸热量急剧减少,造成传热恶化以及堆芯温度升高的现象,这不利于反应堆安全。由于反应堆和冷却剂温度升高,反应堆会进行快速响应,引入负反应性来降低反应堆功率,在这个过程中冷却剂焓和反应堆功率会回落并发生大幅波动,直到恢复下一个稳定状态。

焓的变化曲线证明了流量波动在一定范围内,焓波动量较小,能够保证反应堆的安全,但流量短时间快速变化,会造成焓突变和剧烈波动,对反应堆的安全运行带来挑战,不利于安全运行,因此在实际工程中应避免冷却剂流量急剧变化。

3.3 瞬态工况下压降变化分析

图13为在设置的瞬态工况下,子通道分析程序BRESA的压力变化。观察可知,冷却剂流过燃料组件后的压降变化与流量变化趋势相似,流量增加,燃料组件出口压力减少,即冷却剂在燃料组件中压降增大,当流量从5.553 kg/s增加到11.106 kg/s,压降从11 970 Pa变化为41 315 Pa。在流量突然降低的瞬态过程中,燃料组件压降存在波动,其原因与瞬态工况下焓变化原因近似,冷却剂吸热量的大幅变化,使冷却剂密度也会发生变化,从而导致燃料组件压降也会出现波动。分析可知,冷却剂流量变化对压降影响较大,这会导致燃料组件承受的应力变化较大,因此在流量变化的瞬态工况下对燃料组件的材料要求较高。

图13 燃料组件压力变化

流量变化瞬态计算结果表明,应用子通道方法的S-CO2燃料组件中,冷却剂流量的变化会对压降和焓产生一定的影响,当流量随时间线性缓慢变化时,焓和压降的变化也相对平缓,但当流量在短时间内快速变化时,焓和压降都发生剧烈变化,这对S-CO2反应堆的安全运行带来挑战,发生包壳融毁,燃料组件应力腐蚀等事故的概率增大。为了完善反应堆的安全设计,需要设计更科学高效的容积控制系统,保证冷却剂流量稳定、准确、快速、平稳的进行调节。后续研究将对S-CO2反应堆容积控制系统进行设计和分析,本文的研究内容为S-CO2反应堆系统瞬态运行特性分析提供参考。

4 结论

本文基于Modelica语言建立了S-CO2反应堆子通道分析程序BRESA,对燃料组件内各子通道的流量、密度、焓,以及子通道在轴向上的焓升和温升进行了参数提取和数据分析,得到了稳态运行工况下燃料组件的关键参数,证明了开发的子通道分析程序的合理性,并进行了瞬态工况分析,主要结论如下。

1) 在燃料组件中,采用子通道分析方法,由于子通道模型的横向流动的影响,冷却剂流量会有0.2%～1%的变化量,燃料组件热工参数呈现中间通道焓升较高、边缘通道焓升较低的分布特点,验证了程序的正确性和精确性。

2) 在流量变化的瞬态工况下,当流量线性缓慢降低和升高时,反应堆功率与流量变化趋势相同,但呈现波动变化;焓和出口压力与流量变化呈现相反的趋势,并波动变化。这证明了反应堆在冷却剂流量小幅变化的工况,燃料组件参数在控制系统和流量变化的双重作用下会发生小幅波动变化。该工况下反应堆运行参数的波动在一定范围内,反应堆处于安全状态。

3) 在流量变化的瞬态工况下,流量突然变化会导致功率、焓和压降产生剧烈波动,功率发生剧烈波动并快速下降,焓先急剧升高后大幅波动下降,瞬时焓升达到90～120 kJ/kg,流量突变的瞬态工况对S-CO2反应堆的安全运行带来挑战,不利于反应堆的安全运行。