多孔介质界面对微流道散热器流动与换热性能的影响

陶李春，孙 科

(武汉科技大学省部共建耐火材料与冶金国家重点实验室，湖北 武汉，430081)

近年来，微电子[1]、微燃料电池[2]、微动力[3]和微生物芯片等微尺度科技发展迅猛，任务负荷的增加与有限空间的矛盾推动了电子芯片向小型化、高集成化发展[4]。基于“摩尔定律”预测，电子芯片上的晶体管数目大约每18个月增加1倍[5]，这将大幅提升芯片核心的综合性能。与此同时，芯片内核集成电路密度的增加也导致产热量更大，给电子芯片正常运行的热负荷管理带来严峻挑战。

集成电路产生的高热流密度会使芯片表面温度迅速升高，直接影响芯片的可靠性和寿命，据统计，大约55%的集成电路故障是由温度过高引起的[6]。根据10 ℃法则，接近室温的环境温度每提高10 ℃，电子元件使用寿命将缩短一半[7]。因此，需要对电子发热组件、整体或系统采用合适的冷却技术和结构设计，确保其在所处的环境工况条件下能够可靠、安全、高效地工作。

目前，电子器件的主要散热方式有对流(自然对流、强制对流)、相变换热、射流冲击、喷雾冷却和微流道散热[8]等，其中微流道散热器具有结构紧凑、散热能力强、可批量生产等优点，是高热流密度冷却技术领域的首选。优化散热器材料或微流道结构是提升微流道散热器性能的主要措施，在优化散热器材料方面，通常采用散热性能优异的材料作为微流道散热器的主体，利用涂层等对流固交界面进行表面改性或使用纳米流体等代替传统冷却介质，相关的研究如Calame等[9]制备的碳化硅复合材料微流道较碳化硅微流道温度降低了约30%；
Heidarian等[10]研究了纳米流体在(超)疏水壁面结构微流道内流动的特性，结果表明，在一定条件下采用(超)疏水壁面可有效降低泵送功率；
Chamkha等[11]系统总结了纳米流体在微流道散热器中的应用现状。不过，散热器材料制备难度大、成本高且性能不稳定，导致其发展受到一定限制，通过对微流道结构进行优化设计，同样能达到改善散热器性能的目的。微流道结构优化设计主要有改变横截面积、改变通道形式或加入微结构扰流元等，相关的研究如Gunnasegaran等[12]将横截面为矩形、梯形和三角形的微流道进行对比，总结出换热系数和水力直径间的内在关系；
Zhu等[13]将不同振幅、波长、长宽比的波浪形通道与传统矩形长直通道进行比较，发现换热强化对换热性能的影响远高于压降损失；
Chai等[14]在微流道壁面上引入凸起、凹陷或针肋阵列，有效提高了换热效率及温度均匀性。在针对微流道结构的研究中，大多通过增强扰动的方式来加强换热，但复杂的壁面结构会阻碍冷却工质平滑的流动，导致压降损失升高。有鉴于此，本文建立了散热器壁面结构流动与换热数值模型，在模型正确性得到验证的基础上利用数值计算对不同类型壁面结构的流动与换热性能进行了研究，以期为微流道散热器流动与换热性能的改善提供参考。

图1所示为微流道散热器示意图，散热器整体结构见图1(a)，为了简化计算，以单个通道作为计算域。微流道散热器壁面结构依次为固体壁面、组合壁面和多孔壁面(图1(b))，在该三维计算模型中，计算域的模型外部尺寸L、W、H分别为5000、600、1400 μm，微流道内部高度Hc为1000 μm、宽度Wc为400 μm、水平壁面总厚度Th为200 μm、垂直壁面总厚度Tv为100 μm，多孔介质厚度为Tp(图1(c))。定义α为多孔介质厚度Tp在垂直壁面总厚度Tv中的占比，即α=Tp/Tv。不同类型壁面结构的具体尺寸如表1所示。

表1 壁面结构几何参数

2.1 控制方程和边界条件

在计算域内，假定流体流动是不可压缩的、定常的层流流动[15]，多孔介质是均匀、各向同性的，并且内部流体完全饱和，多孔介质与流体处于局部热平衡状态，忽略重力效应和辐射换热。Klank等[16]指出，Navier-Stockes方程仍然适用于水力直径在1 μm到1 mm之间的微流道。采用Brinkman-Forchheimer扩展达西模型和局部热平衡模型来描述多孔介质内的流动和传热[17]，模型控制方程中的质量守恒方程为：

(1)

动量守恒方程为：

(2)

能量守恒方程为：

(3)

μe=μf/ε

(4)

ke=εkf+(1-ε)ks

(5)

式中，kf为流体导热系数；
ks为固体导热系数。

设置模型边界条件，其中入口处条件为：

(6)

式中，uin是均匀入口速度；
Tin是入口温度。

出口处条件为：

(7)

式中，Pout是出口压力。

下底面处条件为：

(8)

式中，qw是散热器底部热流密度，取值为106W/m2。

固体与液体交界面处条件为：

(9)

式中，Tf是交界面流体温度；
Ts为交界面固体温度；
n是固体壁面的法向量。

其它外部壁面处条件：

-ksTs|n=0

(10)

结合边界条件对控制方程进行求解，获得稳态下计算域内的速度场、压力场和温度场。通过参数化处理可以得到表征微流道散热器流动及换热性能的参数。通常利用摩擦系数来评价散热器的流动性能，摩擦系数f的表达式为：

(11)

式中，Dh是散热器矩形通道的水力直径，Dh=4HcWc/2(Hc+Wc);Δp是压降。对应的雷诺数Re可以表示为：

(12)

泵送功率Ω的表达式为：

Ω=QΔp=uinHcWcΔp

(13)

式中，Q是体积流量。

散热器的换热性能可借助平均努塞尔数Num和总热阻RT进行评价，平均努塞尔数Num的定义为：

(14)

式中，hm是表面换热系数，其计算公式为：

(15)

式中，Tw是通道壁面的面积平均温度，Tf是通道内流体的体积平均温度。总热阻RT的定义为：

(16)

式中，Tmax为换热面的最高温度，A为散热器的底面面积。

2.2 网格划分及独立性验证

基于有限体积法对控制方程进行离散化，采用SIMPLE算法进行迭代求解。在收敛准则方面，控制方程式(1)～式(3)连续两次迭代的相对误差小于10-4即可认为收敛。为了提高数值求解的速度和准确性，采用六面体单元的结构化网格，并在边界层及交界面处对网格进行了细化。在数值求解前，先进行网格独立性验证，采用4种不同数量的网格，对比网格尺寸对数值求解结果的影响，以摩擦系数f作为对比参数，相关结果如表2所示。由表2可见，网格数为432 000时计算所得f值与网格数为960 000时的相应值偏差约为0.5%，后者与网格数为1 228 500时的相应值偏差约为0.07%，均在误差允许范围之内，因此取网格数为960 000进行数值求解。

表2 网格独立性验证

2.3 模型可靠性验证

借助Lee等[18]的实验值对基于本研究所提模型的相应数值计算结果进行验证，实验值与模型数值计算值的对比如图2所示。由图2可知，在相同尺寸和工况下，数值计算结果与实验值基本吻合，误差均在15%以内，因此利用该数值计算模型可以有效预测微流道散热器中的流动与换热情况。

图2 模型可靠性验证

3.1 流动性能

在雷诺数Re为56.8～1534.4的工况条件下，对不同壁面结构的微流道散热器的流动与换热进行数值模拟，散热器的压降损失Δp随雷诺数Re的变化如图3所示。由图3可见，在数值模拟的Re范围内，5组不同壁面结构的散热器Δp都随Re的增加而增大。在同一Re下，α越大则Δp越小。不同壁面结构散热器Δp之间的差异在Re较大时更加明显，例如，当Re为852.5时，组合壁面及多孔壁面散热器的Δp较α为0的固体壁面散热器相应值降低了26.5%～31.4%。Δp的降低不仅意味着摩擦系数f的减小，还有利于减少泵送功率。

图3 雷诺数对压降损失的影响

为了明确组合壁面和多孔壁面改善微流道散热器流动性能的机理，本研究重点考察了壁面中多孔介质对流体速度分布的影响。当Re为852.5时，在5组不同壁面结构的微流道散热器中x为2500 μm处，左侧壁面内部及部分通道局部的速度分布如图4所示。从图4中可以看出，虽然5组微流道散热器几何尺寸相同，但在z方向即垂直于入流方向的截面上，速度分布明显不同。以α为1的多孔壁面散热器为例，约占入口流量2.8%的少量流体从通道渗入两侧的多孔介质中，使通道内的流体速度降低，结合式(11)可知，渗流效应对压降损失降低的贡献较小。另外，根据数值模拟结果，对于α为0的固体壁面，流体在z分别为100、500 μm的交界面处于无滑动边界条件，边界层速度梯度较大，导致流体的黏性切应力大、压降损失高，而对于组合壁面和多孔壁面而言，在相同位置的界面处流速不为0，这表明流体在多孔介质壁面上是滑移的而不是黏在壁面上，这种流动效应与使用超疏水表面进行减阻类似[10]。滑移效应导致交界面处速度梯度较小，压降损失降低。因此，组合壁面和多孔壁面散热器的压降较低是渗流效应和滑移效应共同作用的结果。在组合壁面和多孔介质壁面内部，较高的流速只分布在靠近通道的薄层，如图4中矩形虚线框所示，该薄层厚度约占总壁厚的1/2，当壁面中多孔介质的厚度超过该薄层厚度后，继续增加多孔介质的厚度不能对压降的降低产生较大影响，这也解释了为何图3中随着Re的增加，α为0与α为0.5时的压降损失差距较大，而α分别为0.5、0.75及1时的压降损失却相差无几。

图4 微流道散热器速度分布

3.2 换热性能

图5所示为雷诺数Re对微流道散热器总热阻RT的影响。从图5中可以看出，在数值计算设定的雷诺数范围内，5组不同壁面结构的散热器总热阻RT均随着Re的增加而减小。在同一Re下，组合壁面和多孔壁面散热器的RT总是小于α为0的固体壁面散热器相应值，散热器RT从高到低所对应的α值依次为0、1、0.75、0.25、0.5，这表明壁面中的多孔介质可以增强散热器的换热性能；
当Re大于800时，所有散热器的RT逐渐趋于稳定。以Re为852.5的工况为例，组合壁面和多孔壁面散热器的RT较α为0的固体壁面散热器相应值降低了11.0%～17.4%，这是因为组合壁面和多孔壁面中存在多孔介质，流体可以渗入通道两侧的壁面中，加强了流体与固体之间的换热，分布于散热器底面的热量得到有效消散，换热面的最高温度Tmax降低，根据式(16)，总热阻RT将随Tmax的降低而减小。虽然上述数值计算结果显示固体域为某一特定材料时，5组壁面以α为0.5的组合壁面散热器总热阻为最低，但在使用不同材料作为固体域及对α的间隔进一步细化时，发现散热器总热阻RT最低值所对应的α值会在0.3～0.5之间波动。

图5 雷诺数对总热阻的影响

为了进一步探究散热器总热阻RT与多孔介质厚度占比α之间的内在关联，需考察壁面内部温度的变化趋势。图6所示是雷诺数Re为852.5时，在5组不同壁面结构的散热器中x为2500 μm处沿热流密度方向的壁面内部温度变化趋势。由图6可见，在下壁面处，α为0的固体壁面和α为1的多孔壁面内部温度均高于其它3种组合壁面相应值，其中α为0.5的组合壁面内部温度最低；
在左右壁面处，由于流体渗入壁面加强换热，组合壁面和多孔壁面的内部温度开始降低，其中α为1的多孔壁面因渗入的流体最多，故其内部温度下降幅度最大。不过，这种渗流效应也会带来负面影响，随着多孔介质厚度增加，过多的冷却工质将渗入壁面，根据式(5)可知，左右壁面内部整体的有效导热系数会因此减小，从而在组合壁面和多孔壁面中形成较高的温度梯度，一方面不利于散热器整体的温度均匀性，另一方面上壁面及其附近壁面温度较低会导致该区域与流体换热不足。此外，壁面与流体的对流换热效果也会影响下壁面中的温度分布。综上所述，散热器壁面内部温度分布是渗流效应正面与负面影响达到稳态平衡的结果，单纯的增加壁面中多孔介质厚度占比α，虽然可以降低压降损失，但不利于降低总热阻和提高温度均匀性。图7所示为微流道散热器中雷诺数Re对平均努塞尔数Num的影响。由图7可以看出，对于5组不同壁面结构的散热器，Num均随着Re的增大而增大；
当Re一定时，组合壁面散热器和多孔壁面散热器的Num均高于α为0的固体壁面散热器相应值，这表明壁面中的多孔介质可以加强散热器内部的对流换热，且在Re较大时，加强效果更加显著，在5组不同壁面类型的散热器中，以α为1的多孔介质壁面散热器的Num为最大；
在Re为852.5的工况下，α为0.5的组合壁面散热器和α为1的多孔壁面散热器的Num分别较α为0固体壁面散热器相应值提高了67.8%和86.4%，结合图(6)和式(15)可知，在换热达到稳定状态时，通道下壁面附近流固温差较大，对流换热更强，虽然α为1多孔壁面散热器的Num更大，但是对流换热主要发生在通道下壁面附近区域，这将导致散热器局部温度过高进而影响整体的温度均匀性。

图6 温度分布

在Re为852.5的工况下，固体壁面(α为0)、多孔壁面(α为1)及RT最低的组合壁面(α为0.5)散热器的固体域温度云图如图8所示。由图8可见，固体壁面散热器整体温度明显高于组合壁面散热器及多孔壁面散热器；
多孔壁面散热器中的最大温差为15.9 ℃，其通道入口上部区域的温度较低，通道出口下部区域的温度较高，整体温度均匀性较差，组合壁面散热器中的最大温差为13.9 ℃，整体温度均匀性较高。

图8 温度云图

本文首先构建了微流道散热器壁面结构流动与换热数值模型，在模型的正确性得到验证的基础上利用该模型研究了散热器壁面结构对流动与换热性能的影响。结果表明，当雷诺数Re为852.5时，与多孔介质厚度占比α为0的固体壁面相比，组合壁面和多孔壁面将微流道散热器的压降损失降低了26.5%～31.4%，且α越大，压降损失就越小，这是组合壁面和多孔壁面中因多孔介质存在而产生的渗流效应和滑移效应所致。当Re大于800时，不同壁面结构的散热器总热阻趋于稳定。当Re为852.5时，组合壁面和多孔壁面散热器总热阻较α为0的固体壁面散热器相应值降低了11.0%～17.4%，其中以α为0.5的组合壁面散热器总热阻为最低，并且，与α为0的固体壁面散热器的平均努塞尔数Num相比，α为0.5的组合壁面散热器和α为1的多孔壁面散热器的相应值分别提高了67.8%和86.4%，虽然α为1的多孔壁面散热器的Num较α为0.5的组合壁面散热器相应值稍高，但由于壁面中渗入过多的流体，导致壁面内部有效导热系数降低，局部温度较高，散热器整体温度均匀性较差。在本研究中，壁面中多孔介质厚度占比α为0.5时，散热器的热阻及局部温度最低、平均努塞尔数较高且温度均匀性较好，综合热性能最优，不过，通过使用不同材料作为固体域及对α的间隔进一步细化的结果表明，对于不同的固体域材料，最佳α值会在0.3～0.5之间波动。