印刷电路板换热器中S-CO2换热特性数值分析

2021-02-10刘凯明杨胡朝营赵富龙余霖田瑞峰

哈尔滨工程大学学报 2021年12期

刘凯，明杨，胡朝营，赵富龙，余霖，田瑞峰

(1.哈尔滨工程大学核科学与技术学院，黑龙江哈尔滨 150001； 2.黑龙江省核动力装置性能与设备重点实验室，黑龙江哈尔滨 150001； 3.中国核动力研究设计院，四川成都 610213)

随着核科学与技术的不断发展，核电所占发电总份额越来越多。其中超临界二氧化碳(S-CO2)布雷顿循环系统作为一种新型高效热能转换系统，在气冷堆中占比越来越大。核反应堆内换热器的换热能力对布雷顿循环系统的热效率有着明显的影响，而印刷电路板式换热器(printed circuit heat exchanger，PCHE)作为一种新型微通道紧凑型换热器在高温高压、高热流条件下相比于普通换热器换热性能更加显著[1]。除在流动换热方面，超临界二氧化碳在萃取、喷涂、清洗、药物制备等领域已得到广泛研究、发展和应用[2-3]。二氧化碳在临界区附近物性会发生剧烈变化，尤其是其比热、热导率、粘度和密度等物性在跨临界区域有极大梯度的变化，工质物性对系统运行压力、温度较为敏感，给超临界二氧化碳换热研究带来了很大困难和挑战[4]。学者已经对超临界二氧化碳印刷电路板换热器的换热特性做了许多的实验和模拟计算。徐向阳等[5]通过数值模拟方法对超临界工况下PCHE的传热和阻力性能进行了研究，探讨了数值模型在两侧变物性流体共轭换热工况下的准确性，同时重点研究了不同工况下流体物性导致的压力和温差变化差异。Marchionni等[6]提出了一种用于超临界二氧化碳电力系统中的PCHE建模方法。可以将PCHE模型嵌入完整SCO2功率模块的模型中，以进行优化、瞬态仿真和控制。结果表明，2种建模方法的预测结果吻合良好，证实了一维方法可以用于PCHE的快速仿真和分析。徐哲等[7]采用数值模拟方法计算了跨拟临界点超临界二氧化碳在印刷电路板式换热器内的流动换热特性。结果表明SCO2温度达到拟临界温度时，流通截面内流体温度分布最均匀，因为此时流体的有效导热系数最大；SCO2侧对流换热热阻在总热阻中占比最大。相梦如等[8]采用数值模拟方法对PCHE通道内SCO2-H2O的耦合特性进行了分析，得出了二氧化碳h、Nu随流量的变化趋势。Liu等[9]通过实验分析了直通道PCHE内SCO2-H2O的换热特性，实验结果表明当SCO2从类气体区冷却到拟临界区时，对流传热增强。Chu等[10]通过SCO2-H2O实验平台探究了不同SCO2工作压力下PCHE的传热率和压力损失，结果表明PCHE在较高压力条件下运行时具有更好的综合性能。Chen等[11]以氦气作为工质，通过实验对直通道PCHE的动态响应进行了模拟和分析，结果表明该响应受入口温度变化、氦气质量流量变化以及两者组合的影响。李玮哲[12]分别研究了超临界甲烷在PCHE中的流动换热特性；曾豪等[13]、刘贵军等[14]分别研究了不同通道形状的PCHE换热性能。一些学者也提出了PCHE流动特性的经验关系式。李雪[15]采用计算流体动力学方法对多种角度的Z型通道在Re≤2 000范围内流体的传热与流动规律进行了研究，建立了综合考虑Re、Pr和Zigzag通道角度3因素的传热与流动准则式。

虽然近年来对印刷电路板式换热器的研究不断深入，但对以两侧均为超临界二氧化碳为工质的研究仍不够充分，且研究都集中在较低的雷诺数范围内。本文研究的直通道印刷电路板式换热器冷热两侧均以超临界二氧化碳为工质，在高雷诺数范围内采用数值计算的方法探究不同工况下印刷电路板式换热器换热特性的变化。

1 S-CO2传热模型及计算方法

1.1 几何模型

本文采用平直通道的印刷电路板式换热器作为研究对象，冷热通道分层布置逆流换热，流动工质皆为超临界二氧化碳，如图1所示。PCHE实际结构中还需要在两侧通道连接接口和支架来导入导出工质，在本模型中已经省略，只考虑了芯体部分。该换热器由100个板片组成，每个板片有80个通道；单个管道直径2 mm，板片厚度1.5 mm，管道总长度740 mm。

图1 印刷电路板式换热器基本结构Fig.1 Basic form of printed circuit heat exchanger

本文以一个换热单元为研究对象，上部为热侧通道，下部为冷侧通道，上下左右壁面为周期性边界条件，模型尺寸为3 mm×3 mm×740 mm。

1.2 数学物理模型

文献[16-17]针对微通道内超临界流体的热工水力特性研究成果表明SSTk-ω模型强化壁面处理方法较好，能够对超临界流体给出更准确的流动与换热结果预测，更适用于解决湍流问题。同时，SSTk-ω模型结合了k-ε、标准k-ω模型的优点，即k-ε模型能够对远离壁面的区域给出精准的预测，而标准k-ω模型对近壁面处给出精确的预测，以及模型能够对反压力梯度下流体分离做出高度精确的预测。因此，本文选择SSTk-ω模型进行直通道PCHE流动与换热研究。

CO2物性是基于美国国家标准与技术研究所提出的制冷剂及其混合物热力学和传输特性数据库的材料数据库获得。CO2在拟临界区内物性变化剧烈，尤其是密度、比热容、导热系数和粘度系数靠近拟临界区时有剧烈的陡增或陡降。超临界CO2在类气区和类液区分别符合常规流体性质，在类液区运动粘度及导热系数与温度成反比，类气区运动粘度及导热系数与温度成正比。

由于两侧通道压力不同，CO2物性需要分开设置。考虑到本研究中CO2在PCHE最大压降不超过0.1 MPa，可以将CO2物性作为温度的函数通过分段线性方式输入到FLUENT中，物性分辨率见表1。为了进一步证明分段线性方法得到准确性，将8.5±0.1 MPa和15.5±0.1 MPa在工作温度范围内的物性作对比发现，最大误差不超过2.19%，平均误差在0.25%以内，即分段物性方法合理可行。

表1 物性插入点分布表Table 1 Physical property insertion point distribution

1.3 网格划分及无关性验证

选用ANSYS ICEM软件对三维直通道模型进行结构化网格划分。为了提高计算精度以及充分捕获边界层流动信息，需要对直通道边界层网格进行加密处理。因两侧通道边界层最大不超过0.01 mm所以在划分边界层网格时将两侧通道边界层设置为10层网格，第1层网格高度为0.001 mm，增长因子为1.1。

采用ANSYS FLUENT对模型进行计算，换热器固体材料设置为不锈钢，考虑重力影响，运行压力设置为8.0 MPa以保证CO2处于超临界状态。计算过程中冷热两通道CO2雷诺数Re>4 000，呈湍流状态。

网格数量会对数值模拟结果产生很大的影响，因此为了节约计算资源以及提高计算精度与效率，需要对几何模型进行网格无关性验证。针对构建的几何模型，共划分4种不同数量的网格，分别为480 700、1 019 900、2 832 500和5 491 300。在验证过程中，保持所有模拟工况的边界条件均相同，分别对不同数量的网格进行模拟计算，得到冷热通道压降的变化情况，如图2所示。可以看到，网格数量大于2 832 500时，热通道出口温度和通道压降基本保持不变。因此为了节约计算资源和提高计算效率，选择网格数量为2 832 500的网格模型进行模拟计算与分析。

图2 网格无关性验证Fig.2 Grid independence verification

1.4 模型验证

为了验证数值模型计算结果的可靠性，选用秦娜[18]等计算水平管内SCO2的流动换热实验结果进行验证。其实验装置为长400 mm的正六棱柱形水平管道，通过管外冷却水逆流与交叉流动换热研究了入口压力、质量流量和热流密度等参数对SCO2换热的影响，本文对其实验进行了建模。进行验证计算时入口温度、管道表面热流密度及出口压力的选取与实验条件保持一致。通过计算得到管道壁面热流密度与实验结果进行对比，如图所示。

计算结果如图3所示，数值模拟结果与实验结果趋势一致吻合较好，最大误差低于10%，证明了本文采用的模型及计算方法是可信的。

图3 实验与模拟结果对比Fig.3 Comparison of experimental and simulation results

数值计算与实验存在一定误差的原因：1)由于fluent中CO2的物性采用分段线性方法，并不能完全模拟真实物性；2)实验条件下的热流密度并不是绝对稳定不变的，而是存在小幅度偏离。

2 工况变化对PCHE换热器换热影响

2.1 温度分布

分别改变冷热通道入口温度进行计算，热侧出口压力8.5 MPa，冷侧出口压力15.5 MPa，两侧入口流量均为0.001 4 kg/s，其他参数见表2。

表2 温度变化工况的计算参数Table 2 Calculation parameters for temperature change conditions K

图4为工况5条件下PCHE热侧通道沿轴向的温度云图，由于入口段效应，局部传热系数比其他部分的要高(图5所示)，所以在最靠近入口处截面的中心主流温度范围较广，从壁面到中心的温度梯度较陡；随着流动发展换热进行，截面主流中心温度范围变小，从壁面到中心的温度梯度变缓。

图4 热侧通道沿轴向温度变化云图Fig.4 Contour of temperature change along the axis of the hot side channel

图5 热侧通道换热系数沿轴向的变化Fig.5 Variation of the heat transfer coefficient along the axial direction of hot side channel

2.2 工质温度变化对PCHE换热器换热影响

2.2.1 热侧通道入口工质温度变化

以工况5为基础保持冷侧入口条件不变，改变热侧通道入口温度进行计算得到对流换热系数h和Nu随热侧入口温度变化的曲线，如图6、7。由图可以看出，随着热侧通道入口工质温度的升高，冷热两侧管道换热系数和Nu呈现下降趋势。

图6 对流换热系数随热侧入口温度变化曲线Fig.6 Variation curves of heat transfer coefficient with inlet temperature of hot side

图7 努塞尔数随热侧入口温度变化曲线Fig.7 Variation curves of Nusselt number with inlet temperature of hot side

由于温度增加使CO2热导率上升，流量一定，使得Re减小，因此Nu呈现降低趋势，冷热两侧温差变大使得PCHE热功率增大。在8.5 MPa和15.5 MPa的压力下，热侧和冷侧的超临界二氧化碳温度几乎全部在低密度类气体区，呈现类气体状态，热流密度增加，壁面与工质温差变大(图8)，对流换热减弱[19]。

图8 SCO2与壁面温差随热流密度变化曲线Fig.8 Variation curves of temperature difference between SCO2 and wall surface with heat flux density

2.2.2 冷侧通道入口工质温度变化

以工况5为基础保持热侧入口条件不变，改变冷侧通道入口温度进行计算得到h和Nu随冷侧入口温度变化的曲线，如图9、10。由图可以看出，随着冷侧通道入口温度的升高，冷侧管道换热系数呈现下降趋势，Nu呈现上升趋势；热侧通道换热系数和努塞尔数都呈现下降趋势。

图9 对流换热系数随冷侧入口温度变化曲线Fig.9 Variation curves of heat transfer coefficient with inlet temperature of cold side

图10 努塞尔数随冷侧入口温度变化曲线Fig.10 Variation curves of Nu with inlet temperature of cold side

在15.5 MPa下整个冷侧通道内二氧化碳动力粘度随入口温度的增加而降低，使得Re增大，且Pr减小，Nu呈现升高趋势；冷热两侧工质都处于超临界类气区，入口温度增加，热流密度减小，壁面与工质温差变大，对流换热减弱；在8.5 MPa压力下，二氧化碳热导率随温度升高而变大，使得热侧Nu呈现下降趋势，且冷热两侧温差降低使得PCHE换热功率减小。

2.3 工作压力变化对PCHE换热器换热影响

2.3.1 热侧通道入口工质压力变化

以工况5为基础保持冷侧入口条件不变，改变热侧通道工作压力进行计算得到h和Nu随热侧压力变化的曲线，如图11、12。由图11、12可以看出，随着热侧通道压力的升高，热侧通道h呈现上升趋势，Nu随压力升高先增大后减小；冷侧通道h、Nu都呈现下降趋势。

图11 对流换热系数随热侧压力变化曲线Fig.11 Variation curves of heat transfer coefficient with pressure on the hot side

图12 努塞尔数随热侧压力变化曲线Fig.12 Variation curves of Nu with pressure on the hot side

随着热侧通道压力升高，热侧通道对流换热系数递增，CO2热导率上升，二者共同作用，使得热侧通道Nu随压力先升高后降低。热侧通道换热增强使冷侧工质温度上升，随着温度的上升，CO2导热系数增加，且对流换热系数随压力升高而减小，因此冷侧通道Nu与对流换热系数呈现了相同的变化趋势。

2.3.2 冷侧通道入口工质压力变化

以工况5为基础保持热侧入口条件不变，改变冷侧通道工作压力进行计算得到h和Nu随冷侧压力变化的曲线，如图13、14。由图可以看出，在热通道一侧Nu随压力增加呈现上升趋势，对流换热系数先升高再降低；在冷通道一侧随着压力的升高对流换热系数呈现上升趋势，Nu呈现下降趋势。

图13 对流换热系数随冷侧压力变化曲线Fig.13 Variation curves of heat transfer coefficient with pressure on the hot side

图14 努赛尔数随冷侧压力变化曲线Fig.14 Variation curves of Nu with pressure on the hot side

随着冷侧通道压力升高，冷侧通道对流换热系数升高，CO2热导率随压力升高而增加但努塞尔数呈现上升趋势。CO2热导率对努塞尔数的影响中占主导地位。冷侧通道换热增强使热侧工质出口温度降低，CO2热导率减小，Nu增加使得h出现了先升高后降低的趋势，即在此过程中对流换热系数和热导率的主导地位发生了改变。

对流换热系数和热导率主导地位的变化说明流体层的导热阻力与对流传热阻力的强弱是变化的，增强对流换热及始终让对流换热系数在计算中占主导地位对增强PCHE换热具有重要意义。

根据典型D-B公式，随着热侧通道压力升高，热侧通道内二氧化碳运动黏度减小，Re增大，Pr也增大，热侧通道努塞尔数应随压力升高而增加，但计算结果正与其相反。此类情况说明经典D-B公式并不再适用于超临界CO2，应当作出修正，即考虑壁面温度对传热的影响，引入壁面温度作为相关参数的定性温度。

将温度、压力变化工况的4组计算结果对比可以发现，由于工作压力不同，部分工质工况从类液区向类气区过渡。超临界二氧化碳在类气区和类液区分别符合常规流体性质，在类液区运动粘度及导热系数与温度成反比，类气区运动粘度及导热系数与温度成正比，使得在PCHE换热通道入口温度变化相同大小时h和Nu呈现出不同的变化趋势。

2.4 工质流速变化对PCHE换热器换热影响

2.4.1 热侧通道入口工质流速变化

以工况5为基础保持冷侧入口条件不变，改变热侧通道入口流量进行计算得到h和Nu随热侧流量变化的曲线，如图15、16所示。由图可以看出，随着热侧通道入口流量的增大，热侧通道h递增，冷侧通道h递减，随着流量进一步增大冷侧通道h减小速率放缓。两侧通道Nu也同h一样呈现出相同的变化趋势。

图15 对流换热系数随热侧流量变化曲线Fig.15 Variation curves of heat transfer coefficient with the flow rate on the hot side

图16 努赛尔数随热侧流量变化曲线Fig.16 Variation curves of Nu with hot side flow

随着热侧入口流量增加，工质雷诺数增速较快占主导地位，使热侧Nu呈现明显上升趋势；热侧工质出口温度随流量的增加而上升，工质热导率升高，导致热侧通道h增加。热侧通道流量增加使得冷侧工质温度上升，工质运动黏度增加，Re减小，Pr占主导地位，使冷侧Nu呈现降低趋势；冷侧工质温度升高，导热系数降低，通道换热系数减小。

2.4.2 热侧通道入口工质流速变化

以工况5为基础保持热侧入口条件不变，改变冷侧通道入口流量进行计算得到h和Nu随冷侧流量变化的曲线，如图17、18所示。由图可以看出，随着冷侧通道入口流量的增加，冷热两侧通道Nu、h均呈现递增趋势，但热侧通道增速较缓。

图17 对流换热系数随冷侧流量变化曲线Fig.17 Variation curves of heat transfer coefficient with the flow rate on the cold side

图18 努塞尔数随冷侧流量变化曲线Fig.18 Variation curves of Nu with cold side flow

随着冷侧入口流量增加，冷侧工质Re增速较快占主导地位，使Nu呈现明显上升趋势；冷侧工质出口温度随流量的增加而降低，工质热导率降低，热侧通道h增加。冷侧通道流量增加使得热侧工质温度降低，工质运动黏度降低，Re增加，Pr升高，使热侧Nu呈现升高趋势；热侧工质温度降低，导热系数减小，通道换热系数增加。

在增加冷侧入口流量条件下，热侧通道的换热系数和Nu并没有像增加热侧入口流量条件下的冷侧通道一样呈现降低趋势，这主要是由于CO2在冷热侧不同工作压力下热导率的变化趋势是不同的。热侧流量增加使冷侧出口工质温度升高，整体运动黏度升高，从而Nu、h减小；冷侧流量增加使热侧出口工质温度降低，整体运动黏度降低，从而Nu、h升高，呈现出2种不同的趋势。