贺晓斌++陈文振++肖红光
摘要：流体在稠密栅元内的流动会产生相干结构，而水作为实际反应堆中的冷却剂，探究其在稠密栅内的流动传热将具有很大意义。本文利用URANS（非稳态雷诺应力）法分析不同Re下水在7棒束稠密栅元内流动传热特性。结果表明：在所取的Re下，栅元相同位置处Re对相干结构的强度影响不大，但Re越大，涡结构的速度梯度越大；由于水流动各向异性强，流场和温度分布分布不具规律性，但Re越大，分布梯度越明显；Re越大，沿壁面壁面剪应力和湍动能的最大值越大，达到最大值的位置越滞后，同时随Re增大
壁面剪应力稳定性越好，湍动能稳定性越差。
关键词：稠密栅元 水 相干结构 Re
堆芯冷却剂的流动传热对反应堆的安全和效率有重要影响。当前，稠密栅元（即棒心距与棒径比P/D<1.1）得到了越来越多的关注。采取这种栅元布置方式能够，降低水铀比，使堆芯中子能谱硬化，增大燃料转换比与功率密度，进而提高堆芯效率[1-3]。相关的研究表明：在稠密栅元内，随着流体的流动交混，在横向上会产生大尺度、准周期性、不同于二次流的涡结构，即相干结构，其横向波动速度可达主流速度的10%以上[4-5]。同时，小棒束（如4、7棒束）栅元在工程应用中，中心通道流动受壁面通道影响明显，对其内的相干结构分析会更有代表性。但目前已有的对稠密栅元相干结构的探究都采用空气为工质。本文将在已有研究基础上，采用Fluent CFD方法，探究不同Re下水7棒束稠密栅元内流动传热的特性。
2 几何模型
本文以水为工质。为使流动能达到充分发展状态，栅元由7根长8米的燃料棒组成，采用全通道加热。为了减少变量不同给结果造成的影响，在棒束排列的外围加入与肖红光7棒束模拟相同的六边形围板[6]。流动方向划分如图1。
文中栅元内棒束采用直径D=140mm，P/D=1.06（P为棒芯距），W/D=1.03（W为壁距），棒间距8.4mm，棒壁间距2.1mm，棒束截面如图2所示。
3 网格划分及边界条件
文中网格生成和计算采用六面体结构化网格划分方法，根据所研究的7棒束栅元几何结构的对称性，采用其六分之一通道（即子通道）作为7棒束栅元网格生成的最小单元。整个栅元的网格由子通道网格通过CFD ICEM的阵列功能实现[7]。同时，子通道网格分布分别沿轴向方向、壁面周向方向和垂直于壁面的径向方向进行设置。计算网格如图3所示。
图3计算网格
文中选用肖红光7棒束模拟中的两组Re[4][6]，流动参数如表1所示。选取速度入口和压力出口作为进出口边界条件，网格数为11270784，壁面第一层网格y+<1，由于温升较小，水的物性参数取定值，采用RSMStress-Omega湍流模型[8]。同时，为缩短瞬态计算时间，模拟中先进行稳态计算，稳定后再进行瞬态计算。算例中时间步长指定为0.002 s。采集数据的时间多于一个周期。
表1流动参数
4模拟结果及分析
4.1 相干结构
是否产生相干结构是本次模拟的关键所在。对比流道内3个不同位置的流动数据，如图4所示，点1、2、3分别对应中心棒-外围棒间隙、外围棒-外围棒间隙以及外围棒-围板间隙，线2、3、4分别表示该位置贯穿整个流道。线5是流道中轴线，其端点为中心棒和外围棒圆心构成的等边三角形的中心。
点1-3在Re=38754和Re=64590两种工况下的相干结构如图5，以无量纲的横向速度与入口速度的比值（w/Ub）表征。可见其相干结构形式呈正余弦曲线。
图5（a）、（b）中点1和点2两种工况下相干波动较为紊乱，速度梯度大。Re=64590工况下栅元内流体主流速度大流量多，相干结构使其交混更强，故其横向波动要明显于Re=38754。点1处二者峰值只有0.4%左右，表明此处流体流动以轴向速度为主。点2处所受加热较点1虽少，但由于围板对流动的扰动作用。二者峰值达到12.5%，与横向波动占主流速度达10%以上的实验结果吻合[8]。图5（c）中点3处两种工况下相干波动呈现明显规律性，虽然点3处流量较少，但二者峰值接近，在10%左右，这是由点3处离围板最近，在主流速度较小情况下，围板对流体的阻隔产生的扰动使流体各向速度增大引起。三点处两种工况下涡结构的振幅相近，表明两者相干结构的强度相近，可见在流体物性（尤其是粘性）为定值条件下，同一位置处Re对相干结构的产生没有太大影响，涡结构强度相近，但Re越大，横向波动速度梯度越大。
4.2时均轴向速度
本文以无量纲的轴向时均速度与主流速度的比值作为时均轴向速度的表征。
（a）（b）（c）（d）
图6中心通道无量纲时均流场对比
由于水流动具有各向异性，本文的中心通道选用截面上以x轴为对称轴，中心棒壁面+30?处的中心通道。图6（a）、（c）和（b）、（d）分别给出Re=38754和Re=64590两种工况下+30?中心通道内流场分布等值线图。由图可知中心通道内流场呈梯度分布，且从棒壁面向内逐渐增大。在边界层内，二者流动分布相近，速度比由0增大到0.7，速度梯度最大，离开边界层后，由于从壁面向内越靠近相干涡结构，流体受扰动使流动交混增强，故呈递增趋势且速度梯度在逐步减小，表明从壁面向内相干结构强度在减弱.Re=64590工况下+30?中心通道内速度梯度小于Re=38754工况，而在-30?中心通道内速度梯度大于Re=38754工况。
（a）（b）
文中选用包含上述中心通道的子通道。图7分别为两种工况下子通道内流场分布等值线图。由于水流动的各向异性，等值线分布并不对称。在外围棒间隙区。在子通道内，两者在中心棒周围的中心通道内和外围棒周围的壁面通道内的时均轴向速度呈递进分布，且壁面通道内“凸起”大于中心通道。点1处Re=64590工况下相干结构变化明显，交混强，故其速度比高于Re=38754。在壁面通道内，两种工况下各速度比值相同且Re=64590工況下速度梯度小于Re=38754，这与壁面通道内但Re=64590工况下受扰动小有关。在角通道内，点3处由于两者相干涡结构强度相近且稳定，故其速度分布相似，且由于围板的扰动产生横向波动，故其轴向速度低于主流速度。endprint
4.3 温度
本文引入无量纲温度参数Θ作为温度归一化的标准[6][10]，其表达式为：
式中， 为加热面时均温度； 为主流体时均温度，但由于是均匀加热，故取加热面中心位置截面处温度表示； 为流道内任意位置处时均温度。Θ越大，表示该位置温度越低。
（a）（b）
两工况下栅元子通道无量纲时均温度等值线图如8所示。由图知在二者在边界层内温度变化最大，但在中心通道壁面处Re=64590工况下等温线更靠近壁面，可见Re越大，流动和传热越强，温度梯度大。与速度分布相似，在任一子通道内均有两个对称的“凸起”，且壁面通道内“凸起”大于中心通道，但与速度场分布不同的是从棒壁面向内温度呈逐步增大趋势。中心通道内由于受热多于壁面通道，故其温度也高于壁面通道。在壁面通道内，Re=64590工况下等温线体现出较好对称性，其各向导热更均匀。表明Re越大，栅元内温度分布更均匀，更具有“层次”，整体对称性更好。在角通道内，由于流量和受热最少，二者温度分布相近，且温度都较低。整个子通道内，由于栅元整体温升小，故其等温线最高为0.9，最低为0.7，与初始温度差别不大。
4.4 壁面剪应力
周向剪应力以棒壁上任意位置点处的剪应力与该处时均剪应力的比值进行表征。
两种工况下的无量纲中心通道壁面剪应力随角度分布如图9所示。由图知随着沿壁面角度的增大，两种工况下壁面剪应力整体规律性较差，在φ=0～7.5?范围内，壁面剪应力均有所下降，且Re=38754工况降低速率较大，在φ=12.5～22.5?内，两者变化趋势相似，但Re=64590工况下剪应力低于Re=38754工况。相较而言，Re=64590工况下剪应力分布稳定性较好，其峰值高于Re=38754工况，这与在粘性恒定情况下，流体的流速大，生成相干结构相对强，流动扰动剧烈有关，但其达到峰值的位置约在φ=30?，滞后于Re=38754工况。此外，相对于真实流动，模拟中采用雷诺应力湍流模型，在弯曲壁面会造成一定失真[9]。
4.5湍动能
中心棒壁面的周向湍动能分布如图10所示，以Fluent中定义的湍动能计算方式进行表征。由于Re=64590工况下流体流速和质量流量大，故其各位置处的湍动能要远高于Re=38754工况。在φ=0～7.5?范围内，两者湍动能呈下降趋势，在φ=15～25?内湍动能均逐步增大，随后又稍有降低。整体而言Re=64590工况下湍动能变化梯度更大，稳定性较差，这与Re=64590工况下点1处相干结构较强且变化更剧烈有关，其峰值出现的位置也要滞后于Re=38754工况。
5结论
本文以水为工质对7棒束稠密栅元流动传热特性进行模拟和分析，得出结论如下：（1）相同位置处Re对相干结构的强度影响不大，但Re越大，涡结构的速度梯度越大，结构越不稳定。相干结构是加热面大小、流量及有无外界扰动等因素综合作用的结果。（2）栅元流场和温度分布具有周向异性，且由于水流动各向异性强，其分布不具规律性。但Re越大，分布梯度越明显，对称性也越好。（3）Re越大，壁面剪应力稳定性越好，但湍动能稳定性越差，两者峰值更高但达到峰值的位置要滞后。
参考文献
[1]于意奇，顾汉洋，杨燕华，等. 稠密栅元棒束内流动行为的 CFD 数值模拟[J]. 核动力工程， 2009， 30（2）： 6-11.
[2]鄢炳火，顾汉洋，于雷. 稠密栅元内湍流流体的流动传热[J]. 原子能科学技术， 2011， 45（12）： 1427-1430.
[3]何斯琪，赵冬建，于意奇，等. 7棒束紧密栅元流体流动传热数值研究[J]. 原子能科学技术， 2012， 46（S）： 199-205.
[4]Xiao H G， Yan B H. Numerical study of the turbulent flow and vortex structures in a new lattice[J]. Nuclear Engineering and Design 2014， 275： 344-351.
[5]于意奇，顧汉洋，程旭，等. 不同棒束结构稠密栅元通道内的湍流CFD研究[J]. 原子能科学技术， 2012， 46（4）： 396-403.
[6]肖红光.异形栅元流动传热特性及相干结构数值分析[D].武汉：海军工程大学， 2015.
[7]唐家鹏.ANSYS FLUENT 16.0超级学习手册[M].北京：人民邮电出版社，2016.
[8]Meyer L. From discovery to recognition of periodic large scale vortices in rod bundles as source of natural mixing between subchannels-A review[J]. Nuclear Engineering and Design， 2010， 240（6）： 1575-1588.
[9]曾和义，郭赟. 湍流模型在堆芯热工水力性能分析中的应用[J]. 原子能科学技术， 2011， 45（1）： 20-24.
[10]胡高杰.非均匀壁面粗糙度棒束间流动传热特征的数值分析[D].武汉：海军工程大学， 2013.endprint
