稠密栅元内湍流流体的流动传热

鄢炳火，顾汉洋，于 雷

（1.海军工程大学 核能科学与工程系，湖北 武汉 430033；2.上海交通大学 核科学与工程学院，上海 200240）

当核反应堆内棒束的P／D（P为中心距，D为棒直径）小于1.1时，这种棒束结构通常被称为稠密栅元。稠密栅元内流道间的横向流动较常规棒束间的横向流动更加明显，且呈准周期性的波动，其横向速度波动振幅可达主流速度的10%以上。这种准周期性的流体波动可在很大程度上增强流体的横向交混和传热能力。有关实验和理论分析结果［1－3］均表明，稠密栅元内的这种横向流动的主要特征是大尺度、准周期性的涡结构运动。这种运动受二次流的影响较小，而主要是由湍流大尺度相干结构决定。由于核反应堆堆芯内的冷却剂流动特性与核反应堆系统的热工水力特性密切相关，因而分析稠密栅元内的湍流流动特性、探索稠密栅元内的横向交混特性对新一代核反应堆的设计和安全运行均有重要理论和实践意义。

目前，国外已有一些学者开始研究稠密栅元内准周期性的涡结构运动规律和流动传热特性，文献［3－4］中均利用非稳态雷诺平均模拟（URANS）对这种流动现象进行了研究。但到目前为止，国内尚无这方面的研究。本文对稠密栅元内的湍流流动和传热特性进行分析，首先利用实验数据对计算结果进行验证，然后分析Re和P／D等参数对稠密栅元内的摩擦阻力系数和传热系数的影响。

1 理论模型

1.1 稠密栅元

由于稠密栅元结构具有的一系列优势，在核反应堆系统中，如重水堆和超临界水冷堆均采用这种结构的堆芯布置。Krauss等［1－2］曾对这种稠密栅元中的流体进行了一系列的实验研究。本文以Krauss实验中的稠密栅元为研究对象，其排列方式如图1所示。

图1 稠密栅元Fig.1 Tight lattice

所有热工水力参数均与实验中的参数相同，即：1）实验Re为38 754；2）热工水力直径为33.5mm；3）P／D＝1.06；4）实验工质为常温常压下的空气，其动力粘度为1.919 4×10－5Pa·s；5）热流密度恒定为0.98kW／m2。

1.2 边界条件

在主流方向上，有两类边界条件可供选择，入口／出口边界条件和周期性边界条件。由于周期性边界条件所需的通道长度和网格都更少一些，也可更好地描述充分发展段的流体，因而本文在主流方向上也采用周期性边界条件。栅元轴向长度为实验中所测波长的4倍（600mm）。本文采用的计算域如图2所示。根据实验中横向流体流动的对称性，在径向设定了B－B和C－C两对周期性边界条件。

图2 边界条件Fig.2 Boundary condition

1.3 数值格式与湍流模型

本文采用CFD软件FLUENT6.2进行计算，数值格式如下：1）压力速度耦合方程采用PISO算法进行求解；2）隐式欧拉格式对瞬态项进行离散；3）其它差分格式均为二阶精度；4）利用壁面函数对近壁区的边界层进行求解，第1层网格大小为20＜y＋＜100，主流区采用结构化网格，如图3所示。本文最后得到的网格数为576 000。

图3 网格划分Fig.3 Mesh generation

与直接数值模拟（DNS）和大涡模拟（LES）相比，非稳态雷诺平均模拟（URANS）所需的网格、计算步长和计算资源均更少。因而Merzari等［3］在计算稠密栅元内的流体时均采用URANS模拟。同样本文也采用URANS模拟方法。由于RSM模型是一各向异性的湍流模型，因而可很好地模拟湍流的Reynolds应力和横向流动。本文的湍流模型为非稳态RSM模型。为保证数值残差足够小，计算过程中的Courant数小于0.2。

2 结果分析

本文首先对实验条件下稠密栅元内的流体流动和传热进行模拟，利用Krauss等［1］的实验数据对计算结果进行验证，然后分析Re和P／D等参数对稠密栅元内的摩擦阻力系数和传热系数的影响。

2.1 实验验证

壁面温度和壁面剪应力的分布如图4所示。图4中的参数分别用平均壁面温度和平均壁面剪应力进行了归一化。从图4可看出，URANS的计算结果与实验数据非常吻合，壁面温度与实验结果之间的误差非常小，壁面剪应力与实验值之间的误差也在5%以内。所以本文所选用的湍流模型、数值格式和网格等都是准确、可信的。

2.2 摩擦阻力系数

图4 壁面温度和壁面剪应力的分布Fig.4 Profiles of wall temperature and wall shear stress

图5示出了摩擦阻力系数λ随Re的变化。从图5可看出，利用Blasius阻力系数公式计算得到的摩擦阻力系数与该稠密栅元内的摩擦阻力系数较接近，二者的差异基本上在10%以内。图6示出了Re为38 754（实验Re）时摩擦阻力系数随P／D的变化。传统的理论模型均认为，如果Re相同，则摩擦阻力系数也相同，但图6表明，在稠密栅元内，即便Re相同，P／D的变化也会给摩擦阻力系数带来显著的影响，而传统的摩擦阻力系数计算公式已无法完整地描述稠密栅元内的摩擦阻力系数变化规律，因而也无法完整地描述稠密栅元内的湍流流动特性。

2.3 Nu

图7示出了P／D＝1.06的稠密栅元内的平均Nu随Re的变化。从图7可看出，利用Dittus－Boelter公式计算得到的Nu与该稠密栅元内的Nu较为吻合，二者的差异基本上在15%以内。

图8示出了Re为38 754时不同P／D的稠密栅元内的平均Nu和最小Nu的变化规律。从图8可看出，虽然Re相同，不同P／D的稠密栅元内的平均Nu之间的差异依然很大，显然Dittus－Boelter公式已不能描述Nu的这种变化趋势，因而也无法完整地反映稠密栅元内的传热特性。随着P／D从1.01增加到1.12，平均Nu先增加后减小。P／D＝1.03的稠密栅元内的平均Nu最大，约为170。这说明P／D＝1.03的稠密栅元内的平均传热能力最好。如果核反应堆堆芯内的棒束的P／D＝1.03，此时栅元所占据的空间体积非常小，而传热能力也达到了最佳，这不但有效地改善了系统的传热状况，还可大幅增加核反应堆的功率，进一步提高了系统的效益和经济性。

图7 平均Nu随Re的变化Fig.7 Variation of average Nu with Re

图8 Nu随P／D的变化Fig.8 Variation of Nu with P／D

从图8还可看出，稠密栅元内的最小Nu随P／D的变化与平均Nu随P／D的变化相似。随着P／D从1.01增加到1.12，最小Nu先增加后减小。P／D＝1.03的稠密栅元内的最小Nu最大，约为92。P／D＝1.12的稠密栅元内的最小Nu约为50，远大于P／D＝1.12的稠密栅元内的最小Nu（Nu＝22）。这说明虽然P／D＝1.01的稠密栅元内的平均Nu较大，约为130，但其最小Nu非常小，这说明栅元内的局部传热能力非常差，如果热流密度过大，极有可能发生局部烧毁的现象。

通过上述分析可发现P／D＝1.03是一临界点，P／D＝1.03的稠密栅元内的流动和传热是最安全，也是最高效的。如果核反应堆的P／D在该临界点附近，此时不但系统的传热特性达到了最佳状态，同时可使核反应堆功率达到最大，能大幅提高运行的效益和经济性。

3 结论

本文对稠密栅元内的湍流流动和传热特性进行了分析，计算结果与实验结果较吻合。Re和P／D均会对稠密栅元内的流动传热产生显著影响，但传统的理论模型无法描述P／D对栅元内的摩擦阻力系数和传热系数的影响，因而无法准确地预测稠密栅元内的流动传热特性。P／D＝1.03是一临界点，这种条件下的稠密栅元内的流动和传热是最安全，也是最高效的。此时核反应堆的功率和系统的传热能力可同时达到最大，在保证系统运行的安全性的同时，还能大幅提高系统的运行效益。

［1］KRAUSS T，MEYER L.Characteristics of turbulent velocity and temperature in a wall channel of a heated rod bundle［J］.Experimental Thermal and Fluid Science，1996，12（1）：75－86.

［2］KRAUSS T，MEYER L.Experimental investigation of turbulent transport of momentum and energy in a heated rod bundle［J］.Nuclear Engineering and Design，1998，180（3）：185－206.

［3］MERZARI E，NINOKATA H.Unsteady Reynolds averaged Navier－Stokes simulation for an accurate prediction of the flow inside tight rod bundles［C］∥The 12th International Topical Meeting on Nuclear Reactor Thermal Hydraulics（NURETH－12）.Pennsylvania，USA：［s.n.］，2007.

［4］CHANG D，TAVOULARIS S.Unsteady numerical simulations of turbulence and coherent structures in axial flow near a narrow gap［J］.ASME Journal of Fluids Engineering，2005，127（3）：458－466.