核动力系统分析中细管簇当量元件阻力

2024-02-02杨从新张文兰

排灌机械工程学报 2024年2期

杨从新,张文兰*

(1. 兰州理工大学能源与动力工程学院,甘肃兰州 730050; 2. 甘肃省流体机械及系统重点实验室,甘肃兰州 730050)

在研究核动力系统不同工况下主管路参数变化的过程中,由于核动力系统的复杂性和特殊性,用其本身做试验的成本极高,因此目前国内外研究者普遍采用一维理论形成的半经验半理论方法对系统动态过程进行计算机仿真[1].尽管一维仿真方法能实现快速计算,但该方法计算精度较低,编程工作量大,所以提出建立三维闭式系统模型对内部动态流动进行计算.三维建模过程中,由于核动力系统蒸汽发生器内部细管簇和反应堆压力容器的内部结构复杂,因此提出在不影响主管路内部流动的前提下,将其内部结构简化为等过流截面的U型弯管,并加入合适的局部阻力元件来当量代替[2-3]细管簇或其他结构的流动阻力,从而实现简化三维模型[4]的目的.

影响管道内局部流动阻力的因素很多,如突变截面、渐变截面和阀门等,其中国内外众多学者针对变截面阻力元件进行了大量研究,如胡娟娟等[5]、YIN等[6]、ABED等[7]、AI等[8]通过试验和模拟,分析了各形状孔板的湍流特性,探究了孔板各结构参数与流体流动阻力的关系,及对圆管内部流动的影响.ROSA等[9]、LEE等[10]研究了多级孔板中孔的位置、孔的几何参数以及相邻孔板的安装位置对压降的影响.尹则高等[11]、郝亚丰等[12]分别对渐缩渐扩、渐缩突扩洞塞的水头损失系数进行了研究,并发现水头损失系数随相对长度和半径收缩比的增加而减小.施红辉等[13]研究了不同类型喷嘴结构内部流体力学特性,JIANG等[14]讨论了从层流到湍流速度下,锥直喷嘴内部的流动特性.

孔板结构简单,但截面变化为突变,对管道下游流动影响较大,难以控制总阻力;而线性渐缩渐扩元件又会使圆管边界层分离现象随渐扩角度的增加而加剧[15].另外,由于需要考虑工艺制造和经济性等因素,目前研究者对相连的渐缩渐扩式阻力元件及圆弧渐缩渐扩式阻力元件的流动性能研究较少.因此文中拟根据流体流动规律,设计一种圆弧渐缩渐扩式阻力元件,并且通过数值计算探究该阻力元件对圆管内部流动的影响,得到阻力元件的ζ与L,d/D之间相互影响的规律,从而拟合并验证该阻力元件的ζ与结构参数之间的函数关系式,最后利用该指数函数的计算和预测实现某蒸汽发生器内部细管簇的简化.

1 当量阻力元件设计及阻力系数计算

图1为核动力系统中蒸汽发生器的内部结构示意图.

图1 核动力系统中蒸汽发生器内部结构示意图Fig.1 Schematic diagram of internal structure of steam generator in nuclear power system

为简化其内部细管簇,设计了一种圆弧渐缩渐扩式阻力元件当量代替原模型流动阻力.由于简化蒸汽发生器中U型弯管直径为680 mm,因此将其安置在相同直径的直圆管内进行研究,图2为安装了阻力元件后圆管的二维模型简图,图中D为圆管内径;B为圆管总长度;L为局部阻力元件长度;d为收缩最小直径,设定当量阻力元件的基本结构参数为L和d/D.为了消除进口对ζ的影响,将圆弧渐缩渐扩式当量阻力元件安装在距离进口3倍阻力元件长度处(C=3L)的位置,同时为使出口流态趋于稳定,给定阻力元件到出口的距离大于5L,因此将圆管长度B设置为15 000 mm.

图2 圆管二维模型示意图Fig.2 Schematic diagram of two-dimensional model of circular tube

根据伯努利方程可知,管道1,2截面局部损失Δh近似表示为

(1)

式中:v1,v2分别为对应截面的平均流速,m/s;p1,p2分别为对应截面的平均压强,Pa;g为重力加速度,m/s2;ρ为流体密度,kg/m3;α1,α2分别为对应截面的动能修正系数,由于流态为均匀流,因而取值为1;z1,z2分别为截面对应位置势能,m.

又因为阻力系数ζ与Δh关系为

(2)

根据z1=z2,v1=v2=v得到ζ为

(3)

2 数值计算

2.1 控制方程、边界条件及求解器设置

对当量阻力元件的ζ进行数值计算.设管内为不可压缩定常流动,使用k-ε二方程湍流模型进行数值计算;采用有限体积法离散N-S方程,并通过SIMPLEC算法对速度和压力进行耦合,计算收敛精度设置为10-5,当满足三大守恒方程,监测量进出口总压差稳定时计算收敛.

连续性方程为

(4)

式中:ui为i向的流速分量,m/s;xi为i向的流体流动距离,m.

动量方程为

(5)

式中:p为压强,Pa;μ为动力黏度,Pa·s;μt为湍流黏性系数;t为时间项.

湍动能方程为

(6)

式中:k为湍动能;ε为湍动能耗散率;Gk为由层流速度梯度而产生的湍流动能;σk为湍动能方程的Prandtl数;μtt为t时间湍流黏性系数.

湍动能耗散率方程为

(7)

(8)

(9)

式中:C1ε,C2ε和Cμ为模型常数项.

边界条件设置:入口设置为水平速度等于6.3 m/s的速度进口,出口设置为自由出流出口,壁面采用无滑移壁面边界条件.由于为高Re流动,且所选用湍流模型为k-ε模型,因此采用应用较广泛的标准壁面函数处理黏性底层.

2.2 网格划分和无关性验证

对安装了当量阻力元件的直管进行六面体结构化网格划分,如图3所示.为了准确计算ζ,在阻力元件处进行局部加密,然后对该网格进行无关性验证.理论上,ζ不会随网格数S的变化而变化,因此将ζ作为网格无关性验证的验证量.取L=612 mm,d/D=0.6的计算工况进行验证,验证结果S为424 669,864 864,1 078 064,1 199 536,1 472 016,1 782 840,2 140 160,2 542 248时,ζ为2.17,1.92,1.89,1.90,1.89,1.89,1.90,1.89.

图3 计算域网格Fig.3 Computational domain grid

在计算过程中,网格的密度和精度会对计算结果产生很大的影响。网格数量少会导致计算结果不准确,网格数量太多又会增加计算时间,因此需要对网格进行无关性验证,确保模拟结果不受到网格数影响,从而使结果更加准确.

3 计算结果与分析

为了进一步确定该当量阻力元件是否便于调节和控制系统总阻力,首先对比分析圆弧渐缩渐扩式阻力元件与传统局部阻力元件对圆管内流场的影响.

3.1 圆弧渐缩渐扩式阻力元件对圆管流场的影响

为了更加清晰地观察管内流场变化,在圆管中安装结构参数为d/D=0.5,L=544 mm的圆弧渐缩渐扩式阻力元件进行分析,图4为直管内水平速度vx的分布.由图可知,流体进入阻力元件之前为均匀流态.随流动的发展,流体在靠近阻力元件上游部位逐渐向管道轴线集中,中心轴向速度开始增加,且在阻力元件与上游管道的过渡管壁处均出现低压旋涡区.当流体进入阻力元件后,流体速度在阻力元件最小截面处达到最大.流体流出阻力元件进入下游管道后,逐渐恢复之前的稳定均匀流态.

图4 直管内水平速度的分布Fig.4 Distribution of horizontal velocity in straight tube

与上游管道和阻力元件过渡处相同,下游管道和阻力元件相连位置处也出现了低速旋涡区,2处产生旋涡的原因相同,都是由于管道截面面积突变,使得上游速度发生变化,管壁处流体在进出阻力元件时,水平速度方向突然改变,流体受到阻力元件壁面的阻碍作用产生回流旋涡区.

为了探究圆弧渐缩渐扩式阻力元件对下游流场的影响程度,将其与突缩突扩式阻力元件和线性渐缩渐扩式阻力元件进行对比验证.分析发现,三者在阻力元件进出口与管道相接部分均出现的低速旋涡区面积具有一定差异.从图4中可以看出,圆弧渐缩渐扩式阻力元件的低速旋涡区面积要远小于线性的渐缩渐扩式阻力元件和突缩突扩式阻力元件,这是由于圆弧渐缩渐扩式阻力元件的壁面形状与流体运动轨迹相似,流体过渡平缓,因此产生的流动分离比较轻微.

图5为3种阻力元件引起管内中心轴向速度变化曲线对比,从图中分析可知,圆弧渐缩渐扩式阻力元件中心轴向速度的变化最小,且对下游影响长度最短.这是由于突缩突扩式和线性渐缩渐扩阻力元件截面面积突然改变,流体没有任何缓冲进入阻力元件,与阻力元件壁面产生较大冲击,因此中心轴向速度剧增,下游恢复区长度较长.

图5 中心轴向速度的变化Fig.5 Changes in intermediate axial velocity

综上可得:相同结构参数下的圆弧渐缩渐扩式阻力元件对管道下游的影响长度更短、旋涡区面积更小.因此相同结构参数的变截面阻力元件中,壁面形状与流体流动轨迹相近的圆弧渐缩渐扩式阻力元件更适合当量代替复杂部件内部流动阻力.

3.2 圆弧渐缩渐扩式阻力元件ζ的影响规律

3.2.1Re对圆弧渐缩渐扩式阻力元件ζ的影响

为了排除流体速度对ζ的影响,改变雷诺数Re的大小,验证Re为4 263 498时,ζ只受结构参数的影响.表1为Re在2 984 449～4 689 848时ζ的变化.

由表1可以看出,Re>3.0×106后,改变Re的大小,ζ几乎不变,因此在对应的给定进口速度6.3 m/s下可以直接探讨结构参数d/D和L对圆弧渐缩渐扩式阻力元件ζ的影响.

3.2.2d/D对圆弧渐缩渐扩式阻力元件ζ的影响

对圆弧渐缩渐扩阻力元件进行设计,由于设计元件目的是减小对下游的影响,便于控制总阻力,因此d/D范围设置为0.5～0.9,L长度在0.8～1.1倍D的范围进行取值,探讨d/D对单个阻力元件ζ的影响,具体的计算工况如表2所示.

表2 计算工况Tab.2 Calculated work conditions

图6为d/D与ζ的关系图.

图6 d/D与ζ的关系Fig.6 Relationship between d/D and ζ

图7 工况1—4下直管内水平速度的分布Fig.7 Distribution of horizontal velocity in straight pipes under working conditions 1 to 4

由图6可知,圆弧渐缩渐扩式阻力元件ζ与d/D呈负相关,即d/D越小,流体经过阻力元件流体之间产生的碰撞和摩擦越多,所消耗的能量越多,因而ζ越大.另外,从图中可以看出,d/D越大,L对ζ影响越小,因而d/D对ζ的影响更为显著.

对比安装了d/D=0.5～0.8,L=544 mm圆弧渐缩渐扩式阻力元件后直管内的水平速度vx分布云图如图8所示.d/D对圆弧渐缩渐扩式阻力元件上游管道的影响长度较短,但对管道下游影响长度比较明显,d/D越小,引起管内流体质点运动越混乱,管内速度梯度变化越大,则下游影响长度越长.另外,当流体进入阻力元件中,由于管径的收缩,导致流体在最小直径管壁处形成高速区域,从图中分析发现,该区域面积会随d/D的减小而增大.

图8 L=544 mm的中心轴向速度Fig.8 Central axial velocity of L=544 mm

由图8所示,当流体即将进入阻力元件时,受管壁约束影响,中心轴向速度开始增加,并在管径收缩最小处达到最大值.分析不同d/D引起管内中间轴向速度的变化发现,d/D越小,对轴向速度产生的影响越大,对管道下游影响长度越大,即下游速度的恢复区域长度也越长.同时可以发现,d/D越小,引起的中心轴向速度变化也越大,即流动状态越紊乱,能量损失越大.

湍动能是判断流动是否稳定的一个重要指标.图9为相同L,不同d/D下的湍动能云图,流体在进入管道时湍动能很小,且分布均匀.随着流动的发展,湍动能开始发生变化:在阻力元件内部,管壁处湍动能明显增大,并在阻力元件后半部分的壁面处达到最大,之后流体进入下游管道,湍动能开始减小,最后在下游3D～5D处恢复稳定.对比湍动能云图9发现,管径d/D越小,引起湍动能的变化越大.另外,随着d/D的减小,湍动能向截面扩大位置集中,其原因是随着d/D增大,阻力元件的扩散度增加,流态出现较大混乱．由此可知,d/D越小,下游流场稳定性越差,恢复段越长.

图9 工况1—4直管内湍动能k的分布Fig.9 Distribution of turbulent kinetic energy k in straight pipe under working conditions 1 to 4

3.2.3L对圆弧渐缩渐扩式阻力元件ζ的影响

取L分别为544,612,646,680,714,748 mm的阻力元件进行分析,探讨L对ζ的影响,如图10所示.

图10 ζ与L的关系Fig.10 Relationship between ζ and L

由图10分析可知,相等d/D的阻力元件,L越大,ζ越小.这是因为上游流体进入相等d/D的阻力元件时,L越大,管径的变化越平缓,流体则更难发生流动分离,流态越稳定.而且当流体进入阻力元件截面扩大段,L越大,流入下游管道的流体越平缓,流体紊乱度恢复得越好,ζ相对越小.

3.3 ζ与结构参数函数关系式的推导和验证

根据前面的计算和对比可知,ζ与L,d/D具有一定的规律性,因此首先对其进行相关性分析.由于3个量都为连续性变量,因而采用皮尔逊分析法对L与ζ,d/D与ζ进行双变量相关性分析.

表3为相关性分析结果,皮尔逊相关性系数的大小代表2个变量间相关性的强弱,显著性系数的大小用来判断两变量之间的相关性是否具有显著性.从分析结论中可以看出,ζ与d/D之间的皮尔逊相关性较强,且呈负相关,2个变量间相关性具有显著性;ζ与L显著性大于0.05,该值说明这2个变量之间的相关性较弱,皮尔逊相关性显示其为负相关.

表3 相关性分析Tab.3 Correlation analysis

为了得出圆弧渐缩渐扩式阻力元件ζ与元件本身结构参数之间的规律,建立该阻力元件ζ与L,d/D之间的函数关系式.

根据表3中相关性分析可知,L对ζ的影响较小,因而将该参数忽略,建立起一个ζ与d/D之间的函数关系式.根据ζ随d/D变化的趋势,对数值模拟数据进行非线性曲线拟合,然后筛选误差较小、置信区间高于95%的2个公式如下:

(10)

(11)

以上2个公式均适应于d/D=0.5～0.9时圆弧渐缩渐扩式阻力元件ζ的预测.

为了证明函数关系式的准确性,对拟合关系式进行验证分析.选取满足公式范围内的6个阻力元件,对其进行数值计算得到ζ,然后将数值计算结果分别与2个函数关系式的预测结果进行分析对比,分别得到相对误差e1,e2.

表4为预测值与数值模拟结果的对比.由表可以看出,在忽略L的情况下,2个函数关系式的预测结果均贴近数值计算结果,其中幂函数关系式的预测误差最大为15.9%,指数函数关系式的相对误差最大为26.3%,且2个误差最大点为同一个点.为了更加清晰地对比,给出3种情况下ζ的对比图,如图11所示.从整体上来看,排除误差最大值点,指数函数关系式较幂函数更贴近数值模拟结果,因而选择指数函数为最终函数关系式.