基于ANSYS的高温下热电单偶的焦耳热分布特性研究

沈书林 涂鸣麟 方洋洋 梁清雲

(1. 三峡大学 电气与新能源学院， 湖北 宜昌 443002； 2. 三峡大学 新能源微电网湖北省协同创新中心， 湖北 宜昌 443002； 3. 长江三峡通航管理局， 湖北 宜昌 443002)

基于塞贝克效应的温差发电技术能将低品位的热能直接转化为电能，具有绿色环保、结构简单、坚固耐用、无运动部件、无噪音和使用寿命长等优点[1-3]．温差(即ΔT)是影响热电器件性能的主要因素，大温差意味着较高的热电输出值．在高温工作条件下热电器件两端的温差往往较大，因此提升热电器件在高温工作条件下的热电特性是温差发电技术的研究重点．

热电器件工作时由于自身内阻的存在将不可避免地产生焦耳热，从而影响其热场特性．国内外学者常通过理论分析和仿真模拟的方法对热电器件工作时的热场特性进行研究，如Lineykin等[4]分析了热电机工作时的热场温度分布，用模拟法研究了其热电传热特性；Chen等[5]采用了同样的方法对温差发电装置传热过程进行分析，得到了单元划分影响温度分布模拟准确度的结论；Al-Nimr等[6]探索不断修正温差发电层两端温度值来分析其稳态温度的分布规律．上述研究和设计注意到了热电器件工作时的焦耳热现象，但在等效计算热、冷端通过的热功率时，认为焦耳热在热、冷端均匀分布，而忽略焦耳热的影响．考虑焦耳热对热电单偶热、冷端温度的影响后，以高温硅锗合金热电单偶为例，通过ANSYS仿真高温工作条件下热电单偶的传热过程，研究高温工作条件下焦耳热分布对热电单偶热场特性的影响．

1 热电单偶工作原理

将富空穴的P型热电材料和富电子的N型热电材料的一端相连，形成一个最基本的PN热电单偶(TEGU-Thermoelectric Generation Unit)[7]．当在两种不同半导体材料组成的热电单偶两端分别施以冷源和热源，热、冷端之间就会产生电势差，若外接负载电阻，回路中会产生电流，温差发电工作原理如图1所示．

图1 温差发电的原理示意图

温差发电器件工作时，热电单偶同时受到电场和温度场的共同作用，其热传导过程可描述为[8]

(1)

假设N型电偶臂(N-leg)和P型电偶臂(P-leg)厚度分别为dn和dp，电阻分别为Rn和Rp，热电单偶热、冷端的温度分别为Th和Tc，则N-leg、P-leg稳态条件下的一维热传导方程[9]

(2)

式中，tn、tp分别为N-leg、P-leg的Thomson系数，An、Ap分别为N-leg、P-leg的截面积，Kn、Kp分别为N-leg、P-leg的热导率，dn、dp分别为N-leg、P-leg的厚度，I为回路电流，则T(x)为N/P-leg关于厚度dn/dp的温度分布函数．

传统的温差发电研究中，通常假设热电单偶工作时产生的内部焦耳热有一半反馈给热端、另一半由冷端导出[10]，这样可以得到热电单偶热端注入热功率Qh和冷端导出热功率Qc的近似结果[11]

(3)

式中，rtot为热电单偶内部总电阻，ΔT为热电单偶两端的温差．

2 基于ANSYSY的热电单偶建模与仿真

2.1 热电单偶建模

热电单偶两电偶臂通常使用电导率较高的金属导体(如Cu片)连接，由于结构和制造工艺的限制，金属导体与电偶臂之间的接触效应会对热电单偶的温度场分布产生较大的影响[12]，因此在ANSYS物理模型中需要考虑接触热阻和接触电阻等影响因素．考虑接触效应后，热电单偶总内阻rtot和总热阻ztot可由下式表示[13]

(4)

式中，ri为P/N-leg内电阻，rc为P/N-leg和金属导体之间的接触电阻，A=An+Ap为热电单偶截面积，ρpc、ρnc分别为P-leg、N-leg与金属导体平均接触电阻率；λp、λn分别为P-leg、N-leg热传导系数，zpc、znc分别为P-leg、N-leg与金属导体等效接触热阻．

硅锗合金是目前研究比较成熟的高温热电材料，在700 K以上有着良好的热电性能，并且在T=1 200 K时，优值ZT=1[14]．以高温热电材料锗化硅为例，N-leg为Si70Ge30P3(30%)和Si95Ge5(70%)的混合物，P-leg为Si70Ge30B5(25%)和Si95Ge5(75%)的混合物[14]，金属导体为Cu片，热电单偶模型的物理尺寸见表1．

表1 元件尺寸

实际中温差热电单偶和Cu片之间并不存在接触层，接触层的作用是代替接触热阻和接触电阻，使计算记过更接近实测值．ANSYS中物理模型如图2所示．

图2 物理模型图

2.2 边界条件与初值

通过ANSYS软件定义N-leg、P-leg、接触层、Cu片及负载电阻的属性参数，热电单偶的材料属性参数见表2[15-16]．其中N-leg、P-leg、接触层、Cu片的单元类型为Thermal Electric、Brick 8node 69，即8节点六面体单元；负载电阻R0为Circuit124单元[17]．

表2 材料属性参数

设置单元格为0.001，通过Mesh操作进行网格划分；设置热端Cu片表面热对流系数为h1=200 W/(m·K)，热源温度为T1=1 200 K，冷端N-leg和P-leg表面热对流系数为h2=1 000 W/(m·K)，冷源温度为T2=300 K，负载电阻R0=4×10-4Ω，工作回路电流根据式(2)设置．

2.3 仿真计算

由于ANSYS自带的温差发电单元Solid226假设一半反馈给热端，另一半由冷端导出[18]，热端注入热功率和冷端导出热功率由式(2)得到，在高温工作条件下误差较大，故采用ANSYS电-热耦合模块，工作回路电流I需要在自定义函数功能Functions中单独设置，ANSYS分析计算流程如图3所示．

图3 ANSYS分析流程图

由于模拟的是热电单偶在高温环境下受到电场和温度场共同作用而产生的电-热耦合效应，从初值到最终的稳定解需要进行多次计算，运用APDL命令循环语句“do while”进行迭代运算．图4为迭代10次后的热电单偶温度场分布云图．

图4 热电单偶温度场分布云图

3 热电单偶中纵向温度分布研究

为研究高温热源对热电单偶温度分布的影响，固定冷源温度T2=300 K，设负载电阻R0=0，改变热源温度T1，稳定下X=0、Z=0上所有节点纵向(Y轴负方向)的温度分布如图5所示．

图5 R0=0时热源温度对热电单偶纵向温度分布的影响

由图5可知，固定冷源温度T2，随着热源温度T1的降低，热电单偶内部焦耳热逐渐减小，热电偶臂内部温度分布线性化愈加明显．

为进一步分析热电单偶内部焦耳热大小和温度分布之间的关系，固定热源温度T1=1 200 K、冷源温度T2=300 K，记负载电阻R0与热电单偶内阻的比值为m，改变m值，则热电单偶纵向温度分布的变化如图6所示．

图6 T1=1 200 K时负载电阻对热电单偶纵向温度分布的影响

根据图6，将m为无穷大(即热电单偶空载)时的温度分布曲线与存在内部焦耳热时进行比较可以得出，由于内部焦耳热的存在，热电单偶热、冷端的温度均升高，但两端的温差值减小；随着m值的增大，热电单偶内焦耳热减小，其内部温度由非线性分布逐渐转为线性分布，热电单偶空载时内部焦耳热Qr=0，其内部温度完全呈线性分布．

由图5和图6，接触层内温降显著，可见接触效应对热电单偶温度分布的影响较为明显，说明模型中接触层的设置是有必要的．

4 高温工作条件对焦耳热分布的影响

根据热力学第二定律，热量总是从温度高的物体传递到较冷的物体[19]，热电单偶实际工作时热端的温度要高于冷端，故在第三类边界条件下热端注入的热功率和内部产生的焦耳热均由冷端导出，但为了简化计算，传统的温差发电研究中在等效计算热电单偶热端注入热功率和冷端导出热功率时，通常假设产生的焦耳热一部分反馈给热端、另一部分由冷端导出．

定义由冷端导出的焦耳热与总焦耳热的比值为焦耳热分布系数s，则反馈给热端的部分占总焦耳热的比值为(1-s)，s的计算式为

(5)

改变热源温度T1，冷端热对流系数h2分别取200 W/(m·K)、400 W/(m·K)和1 000 W/(m·K)，s的拟合曲线如图7所示．

图7 热源温度对焦耳热分配系数s的影响

由图7可知，在T1<500 K的中低温条件下，热电单偶产生的内部焦耳热较小，等效计算中由冷端导出的热功率与反馈给热端的热功率接近，式(4)中关于热端注入热功率Qh和冷端导出热功率Qc的计算误差不大；在T1>900 K的高温工作条件下热电单偶两端的温差较大，内部产生的焦耳热也急剧增大，从冷端导出的热功率要比反馈给热端的焦耳热功率大得多；当T1>1 200 K时，内部产生的焦耳热几乎全部由冷端导出．

随着热源温度的升高，焦耳热分布系数s值逐渐增大，越来越多的焦耳热功率由冷端导出，直至全部由冷端导出．传统温差发电模型中，认为工作时产生的内部焦耳热有一半反馈给热端、另一半由冷端导出的假设已经不再成立．由式(3)计算Qh和Qc，在高温工作条件下势必带来较大的的偏差．

5 结 论

在考虑接触效应的基础上，通过ANSYS软件对高温工作条件下锗化硅合金热电单偶进行仿真，分析了其焦耳热分布特性，获得了稳态温度分布曲线，最后分析了热源温度对焦耳热分布的影响，得出以下结论：

1)空载时，热电单偶内部温度为线性分布；外接负载电阻时，由于焦耳热的影响，热电单偶内部温度由线性分布变为非线性分布．

2)随热源温度的升高，反馈给热端的焦耳热占比不断减小，冷端的焦耳热占比不断增大，直至焦耳热全部由冷端导出．

3)在高温工作条件下，传统温差发电模型中认为反馈给热端的焦耳热和由冷端导出的焦耳热各占一半的假设不再成立．