电臂尺寸和安装条件对热电模块力学性能的影响

Effect of Thermoelectric Leg Size and Mounting Conditions on Mechanical Properties of Thermoelectric Modules

摘要

针对热电制冷模块中热电臂尺寸对模块整体热应力场的影响，以C-31106型热电制冷模块为研究对象，建立由31对热电偶组成的热-电-力学耦合多物理场模型.采用有限元仿真方法，考虑钎料层的存在并引入表征黏塑性材料的Anand修正本构方程，建立3种安装条件下的热电模型，分析不同电流下热电臂高度和截面边长对热电模块温度场和力学性能的影响规律.结果表明，当热电臂处于最佳高度时，在给定范围内热电臂截面边长越大，越有利于提高热电模块力学可靠性，当工作电流为2.6 A，截面边长从0.99 mm增加到1.09 mm时，冷、热端的最大von Mises应力分别减小了74.67%和22.60%.对模块冷端和热端分别采用固定约束和机械压装的安装方式可有效提高热电模块的力学可靠性，当采用上述安装方式时，冷、热端的最大von Mises应力相比两端只采用固定约束的安装方式分别减小了48.76%和74.02%.

Abstract

Aiming at the influence of the thermoelectric leg size on the overall thermal stress field of the module in thermoelectric refrigeration modules， this paper takes the C-31106 thermoelectric refrigeration module as the research object and establishes a thermoelectric-mechanical coupled multiphysics field model composed of 31 pairs of thermocouples. Using the finite element simulation method， considering the existence of brazing layers and introducing the Anand modified constitutive equation characterizing the viscoplastic material， the thermoelectric model under three mounting conditions is established， and the influence laws of the height and the cross-section side length of thermoelectric legs on the temperature field and mechanical properties of the thermoelectric module under different currents are analyzed. The results show that when the thermoelectric leg is at the optimal height， the larger the cross-section side length of the thermoelectric leg is （within a given range）， the more favorable it is to improve the mechanical reliability of the thermoelectric module. When the operating current is 2.6 A and the cross-section side length is increased from 0.99 mm to 1.09 mm， the maximum von Mises stresses at the hot and cold sides are reduced by 74.67% and 22.60%， respectively. The mechanical reliability of the thermoelectric module can be improved by using fixed constraints and mechanical press-fit mounting for the cold and hot sides of the module， respectively. The maximum von Mises stresses on the hot and cold sides of the module are reduced by 48.76% and 74.02%， respectively， compared with the fixed restraints on both sides.

Graphical abstract

关键词

热电制冷 / 制冷器 / 安装条件 / 黏塑性材料 / 热-电-力学耦合 / 有限元法

为了适应航天航空^［1］、航海^［2］、军事^［3］等领域任务需求的提高和任务类型的细化，高度集成化、智能化的电子设备应用日益广泛. 电子系统的高度集成化和微型化，导致器件的组装密度和复杂程度越来越高，系统内部的功耗和发热量剧增，会在局部产生热点，对器件的热稳定性和热可靠性造成不利影响，严重时会使器件直接失效，以上问题对芯片等电子元件的散热提出了更高要求. 传统的热管理方式因安装体积大、维护不便、噪声大等原因^［4］在芯片高效散热上都存在一些缺点或不足，亟须发展占用空间小、维护简单、高效的热管理技术.

热电制冷技术具有体积小、便于集成化和模块化等特点，相比传统的制冷方式，具有控制灵活、无制冷工质、无噪声、无振动、制冷迅速等优点，在精密温度控制^［5］、小空间制冷^［6］、微电子^［7］等领域已得到广泛应用. 如精密温控方面，王镇锐等^［8］利用空间站货运飞船天舟一号作为微重力流体实验平台，开发出一套半导体制冷与泵驱两相回路相结合的精密温控系统，并通过天地数据对比，实现了-5~40 ℃内±0.2 ℃的精确温控. 宋雪冬等^［9］设计了一款采用半导体制冷的高精度星敏感器，通过实验验证，其控温精度可达±0.25 ℃.

在热电材料研究方面，目前的研究热点是提高热电材料的优值系数（Z值），进而提高热电转换效率^［10］. 一些学者^［11-12］在已有热电材料的基础上，通过采用不同制备方法，例如电化学反应沉积法、磁控溅射法、离子烧结法、区熔法、粉末冶金法等来改善材料的热电性能. 近年来随着材料理论和制备工艺的快速发展，热电材料的综合性能已得到显著提升. 但新型热电材料在工程应用时仍然面临制备成本较高、制备工艺复杂和材料热稳定性不足等问题.

在热电装置的热力学分析和优化方面，Chen等^［13］提出了一种折算面积热阻和制冷率密度的分析方法. Baldry等^［14］利用热电制冷器设计了一种高性能的散热器，通过改变结构参数将散热温度降低了11.7 ℃. Zhu等^［15］以熵产最小化为目标建立热力学模型，分析并优化热电制冷器散热性能. Meng等^［16］通过引入有限时间热力学，对结合翅片热沉的热电制冷器模型进行分析和优化. 王子成等^［17］采用有限元仿真和实验相结合的方式，得到了热电制冷器有效工作电流、有效热负荷和有效制冷系数的范围. 孟凡凯等^［18］在绝热表面空间和非绝热表面空间内建立了热电制冷器计算模型，得到了密闭空间热电制冷器工作参数的瞬态变化规律与瞬态特性. 国外的研究发现^［4，19］，减小热电臂高度、增大热电臂截面积会使热电制冷器具有更好的冷却性能，通过减小热电臂的高度，制冷量增加了6.24 W，通过增大热电臂的截面积，制冷系数（Coefficient of Performance，COP）提高了53%.

综合以上研究，有关热电装置的研究主要集中在提高热电转化性能上^［20］，热电装置的力学可靠性研究尚存在不足. 在考虑热、电、力三场作用的同时，对热电模块进行多物理场建模和分析，对于提高热电装置力学可靠性具有重要意义. 现有热电结构-应力文献^［21-24］在分析热电制冷模块力学性能时，只关注了热电臂、电极片、陶瓷片等主要器件的几何参数变化对局部应力分布的影响，缺少对整体应力水平的影响分析. 周林等^［25］从模块内钎料层的角度，分析了钎料层作为黏塑性材料对模块整体力学性能的影响，但没有进一步分析其他内部器件对模块整体力学性能的影响. 在实际应用中^［26］，热电制冷模块与热电发电模块相比，服役工况更为复杂. 复杂的工况要求热电模块不仅制冷性能出色，而且要有较高的力学稳定性^［27-28］，从而保证热电模块具有长期的可靠性. 因此，寻找合适的模块安装方式以及优化内部器件尺寸是提高模块力学可靠性的关键.

近年来，有限元法（Finite Element Method，FEM）在结构力学、热力学和传热学等领域均得到了广泛应用^［29］. 有限元法的优势是通用性强，可建立任意形状的结构模型，适用多种材料和复杂的边界条件. 有限元软件COMSOL可支持热、流体、电磁、结构、声学及多物理场的耦合，既能计算三维立体模型，又能适应复杂的热-电-力学耦合分析，相比其他有限元软件，在多物理场耦合分析方面具有突出优势. 本文利用有限元法，采用COMSOL软件，建立热电模块的热-电-力学多物理场耦合模型，在不同安装条件和电流工况下，对热电臂的高度和截面边长进行变参分析，以期获得热电臂尺寸和安装条件对热电制冷模块整体温度场和热应力场的影响规律，并得到热电制冷模块处于最佳力学可靠性时热电臂的结构尺寸和合适的安装条件.

1 模型多物理场微分方程

1.1 热电有限元方程

在热电制冷计算中，内部热传导的方程可以表示为：

\rho C\frac{\partial T}{\partial t}+\nabla \cdot \mathit{q}=\dot{q}（1）

式中：ρ为物质密度；C为物质比热容；T为热力学温度，\dot{q}为单位体积产热率； q 为热流向量.

电流连续性的方程可以表示为：

\nabla \left(\mathit{J}+\frac{\partial \mathit{D}}{\partial t}\right)=\mathrm{0}（2）

式中： J 为电流密度向量； D 为电荷密度向量.

热-电耦合的连续性方程可以表示为：

\mathit{q}=\mathit{P}\mathit{J}-\mathit{K}\nabla T（3）

式中： P 为珀耳贴系数矩阵；K为导热系数矩阵.

\mathit{J}=\left(\mathit{E}-\mathit{\alpha }\nabla T\right)\mathit{\sigma }（4）

式中： E 为电场强度向量； α 为塞贝克系数矩阵； σ 为电导率矩阵.

电介质的方程可以表示为：

\mathit{D}=\mathit{E}\mathit{\epsilon }（5）

式中： ε 为介电常数矩阵.

当处在无时变磁场环境时，无旋电场中的 E 可以表示为：

\mathit{E}=-\nabla \phi（6）

式中：φ为电势. 将式（3）~式（6）与式（1）和式（2）进行耦合化简后，即可得到热-电耦合的方程.

\mathit{q}=P\cdot \left(\mathit{\sigma }\cdot \left(\nabla U\right)+\mathit{J}\right)（7）

\mathit{P}=\left|\mathit{\alpha }\right|T（8）

\mathit{J}=-\mathit{\sigma }\mathit{\alpha }\nabla T（9）

在不同工作电流下，计算出式（7）~式（9）的结果后即可获得模块的温度场和电场分布.

1.2 热应力有限元方程

由于热电模块所用的碲化铋材料为各向异性材料^［30］，材料中温度梯度的变化会直接影响其导热系数和电导率. 材料的瞬态导热方程可以表示为：

\begin{array}{l}\rho c_{\mathrm{P}}\frac{\partial T}{\partial t}=\frac{\partial }{\partial x}\left[k\left(z\right)\frac{\partial T}{\partial x}\right]+\frac{\partial }{\partial y}\left[k\left(z\right)\frac{\partial T}{\partial y}\right]+\\ \frac{\partial }{\partial z}\left[k\left(z\right)\frac{\partial T}{\partial z}\right]\end{array}（10）

式中：c_P为材料的定压比热容；k（z）为法向方向上的热导.

要对模块进行应力分析，还需要利用柯西（Cauchy）方程的位移-应变关系^［31］，3个方向的位移-应变关系可以分别表示为：

{\overline{\epsilon }}_{xy}=\frac{\mathrm{1}}{\mathrm{2}}\left(\frac{\partial \overline{u}}{\partial \overline{y}}+\frac{\partial \overline{v}}{\partial \overline{x}}\right)（11）

{\overline{\epsilon }}_{yz}=\frac{\mathrm{1}}{\mathrm{2}}\left(\frac{\partial \overline{v}}{\partial \overline{z}}+\frac{\partial \overline{w}}{\partial \overline{y}}\right)

{\overline{\epsilon }}_{zx}=\frac{\mathrm{1}}{\mathrm{2}}\left(\frac{\partial \overline{u}}{\partial \overline{z}}+\frac{\partial \overline{w}}{\partial \overline{x}}\right)

式中：\left(\overline{u},\overline{v},\overline{w}\right)为位移分量，表示变形后新位置在x-y-z轴上的投影；{\overline{\epsilon }}_{ij}为应变张量.

借助非对称矩阵和牛顿方法可得位移和应力间的关系为：

\left\{\begin{array}{c}{\overline{\sigma }}_{xx}\\ {\overline{\sigma }}_{{}_{yy}}\\ {\overline{\sigma }}_{zz}\\ \begin{array}{c}{\overline{\sigma }}_{yz}\\ \begin{array}{c}{\overline{\sigma }}_{z\mathrm{x}}\\ {\overline{\sigma }}_{xy}\end{array}\end{array}\end{array}\right\}=\frac{\overline{\mathit{E}}}{\left(\mathrm{1}-\mathrm{2}\nu \right)\left(\mathrm{1}+\nu \right)}\left[\begin{array}{cccccc}\mathrm{1}-\nu & \nu & \nu & \mathrm{0}& \mathrm{0}& \mathrm{0}\\ \nu & \mathrm{1}-\nu & \nu & \mathrm{0}& \mathrm{0}& \mathrm{0}\\ \nu & \nu & \mathrm{1}-\nu & \mathrm{0}& \mathrm{0}& \mathrm{0}\\ \mathrm{0}& \mathrm{0}& \mathrm{0}& \mathrm{1}-\mathrm{2}\nu & \mathrm{0}& \mathrm{0}\\ \mathrm{0}& \mathrm{0}& \mathrm{0}& \mathrm{0}& \mathrm{1}-\mathrm{2}\nu & \mathrm{0}\\ \mathrm{0}& \mathrm{0}& \mathrm{0}& \mathrm{0}& \mathrm{0}& \mathrm{1}-\mathrm{2}\nu \end{array}\right]\left\{\begin{array}{c}\begin{array}{c}\begin{array}{c}{\overline{\epsilon }}_{xx}\\ {\overline{\epsilon }}_{yy}\end{array}\\ {\overline{\epsilon }}_{zz}\end{array}\\ {\overline{\epsilon }}_{yz}\\ {\overline{\epsilon }}_{zx}\\ {\overline{\epsilon }}_{xy}\end{array}\right\}-\left\{\begin{array}{c}\begin{array}{c}\begin{array}{c}\mathrm{1}\\ \mathrm{1}\end{array}\\ \mathrm{1}\end{array}\\ \mathrm{0}\\ \mathrm{0}\\ \mathrm{0}\end{array}\right\}\frac{\overline{\mathit{\alpha }}\overline{\mathit{E}}\overline{T}}{\mathrm{1}-\mathrm{2}\nu }（12）

式中：ν为泊松比.

热电模块在运行过程中，会产生温度梯度，由于材料热延伸率的原因，材料在各个方向上受压或拉伸，符合三向应力状态下塑性或脆性材料普遍发生屈服失效的情况^［22］. 因此，在分析热电模块的热应力时，采用第四强度理论^［32］. 第四强度理论的强度条件可以表示为：

{\sigma }_{\mathrm{r}}=\frac{\mathrm{1}}{\sqrt[]{\mathrm{2}}}\sqrt[]{{\left({\sigma }_{\mathrm{3}}-{\sigma }_{\mathrm{1}}\right)}^{\mathrm{2}}+{\left({\sigma }_{\mathrm{2}}-{\sigma }_{\mathrm{3}}\right)}^{\mathrm{2}}+{\left({\sigma }_{\mathrm{1}}-{\sigma }_{\mathrm{2}}\right)}^{\mathrm{2}}}\le \left[\sigma \right]（13）

式中：σ_r为冯·米塞斯（von Mises）应力；σ₁、σ₂和σ₃为构件危险点处的3个主应力；［σ］为材料许用应力值.

1.3 Anand修正本构方程

为了更准确地计算热电模块的力学性能，本文将模块内的钎料层考虑在内. 根据已有文献^{［25，33］}，目前被广泛应用于微电子焊接的钎料是无铅锡膏，其中，SAC305是最常用的一种锡膏，为了准确模拟钎料层，可用Anand修正本构方程^［34］来描述钎料层的力学性能. Anand修正本构方程主要用于表征无铅锡钎料层的弹性和非弹性形变，其塑性应变率可用流函数描述为：

{\dot{\epsilon }}_{\mathrm{p}}=A\mathrm{e}\mathrm{x}\mathrm{p}\left(-\frac{Q}{RT}\right){\left[\mathrm{s}\mathrm{i}\mathrm{n}\mathrm{h}\left(\xi \frac{\sigma }{s}\right)\right]}^{\mathrm{1}/m}（14）

式中：{\dot{\epsilon }}_{\mathrm{p}}为非弹性应变率；A为常数；R为通用气体常数；T为热力学温度；Q为活化能；ξ为应力乘子；m为应变灵敏指数；s为方程的内部变量.

内部变量s可通过演化，得到动态应变硬化和恢复的演化方程.

\dot{s}=\left[h_{\mathrm{0}}{\left(\mathrm{1}-\frac{s}{S^{\mathrm{*}}}\right)}^a\mathrm{s}\mathrm{i}\mathrm{g}\mathrm{n}\left(\mathrm{1}-\frac{s}{S^{\mathrm{*}}}\right)\right]{\dot{\epsilon }}_{\mathrm{p}}（15）

S^{\mathrm{*}}=\begin{array}{c}Š\end{array}{\left[\frac{{\dot{\epsilon }}_{\mathrm{p}}}{A}\mathrm{e}\mathrm{x}\mathrm{p}\left(\frac{Q}{RT}\right)\right]}^n（16）

式中：a为硬化/软化相关的应变率敏感系数；h₀为硬化/软化常数；S^* 为s与温度和应变速率有关的饱和值；Š为饱和系数；n为变形抗力饱和值的应变敏感指数. 由上述Anand修正本构方程可知，当确定了A、Q、Š、h₀、S₀、m、n、ξ、a这9个参数后（S₀为变形阻力饱和系数），即可表征黏塑性材料的热力学性能.

2 模型建立

2.1 几何模型

建立由31对热电偶组成的热电制冷模块，如图1所示.其中L为热电臂截面边长，H为模块高度，h为热电臂高度.为使仿真结果与实际情况一致，参考某公司的C-31106型热电制冷模块，长、宽均为15.30 mm，高度为3.70 mm；热电臂的长、宽均为1 mm，高度为1.86 mm. 热电模块内部的热电偶结构如图2所示，一对热电偶主要由p型和n型热电臂、钎料层、铜电极片、上下陶瓷片组成.

2.2 材料设定

要对运行中的热电模块进行静力学分析，须定义钎料层为黏塑性材料，将模块其余部分定义为线弹性材料. 当钎料层采用SAC305时，可采用Anand修正本构方程^［34］来表征其热力学性能，Anand修正本构方程参数如表1所示. 通过元素掺杂的方式调整能带结构，降低热导率，最终实现提高热电材料的优值系数Z值^［35］. 本文所涉及的材料物性参数如表2所示^［36］. 由于p型热电臂和n型热电臂材料的热导率、电阻率是随温度变化的物性参数，因此，p型和n型热电臂的变物性参数拟合公式^［36］分别为：

k_{\mathrm{p}}(T)=\mathrm{1.47}-\mathrm{190}\left(T-\mathrm{300}\right)\times {\mathrm{10}}^{-\mathrm{5}}+

\mathrm{1.99}{\left(T-\mathrm{300}\right)}^{\mathrm{2}}\times {\mathrm{10}}^{-\mathrm{5}}（17）

k_{\mathrm{n}}(T)=\mathrm{1.64}-\mathrm{161}\left(T-\mathrm{300}\right)\times {\mathrm{10}}^{-\mathrm{5}}+

\mathrm{2.56}{\left(T-\mathrm{300}\right)}^{\mathrm{2}}\times {\mathrm{10}}^{-\mathrm{5}}（18）

{\rho }_{\mathrm{p}}(T)=\mathrm{8.8}\times {\mathrm{10}}^{-\mathrm{6}}+\mathrm{51.9}\left(T-\mathrm{300}\right)\times {\mathrm{10}}^{-\mathrm{9}}+

\mathrm{79}{\left(T-\mathrm{300}\right)}^{\mathrm{2}}\times {\mathrm{10}}^{-\mathrm{12}}（19）

{\rho }_{\mathrm{n}}(T)=\mathrm{8.2}\times {\mathrm{10}}^{-\mathrm{6}}+\mathrm{38.7}\left(T-\mathrm{300}\right)\times {\mathrm{10}}^{-\mathrm{9}}+

\mathrm{22}{\left(T-\mathrm{300}\right)}^{\mathrm{2}}\times {\mathrm{10}}^{-\mathrm{12}}（20）

{\alpha }_{\mathrm{p}}(T)=\left(\mathrm{22}\mathrm{224.0}+\mathrm{930.6}T-\mathrm{0.990}\mathrm{5}{T}^{\mathrm{2}}\right)\times {\mathrm{10}}^{-\mathrm{9}}（21）

{\alpha }_{\mathrm{n}}(T)=-\left(\mathrm{22}\mathrm{224.0}+\mathrm{930.6}T-\mathrm{0.990}\mathrm{5}{T}^{\mathrm{2}}\right)\times {\mathrm{10}}^{-\mathrm{9}}（22）

2.3 边界条件和计算方法

为了简化计算，在如图1所示的模型基础上作出模型假设并设置如下边界条件：

1）除冷端和热端外，所有表面都视为绝热表面.

2）由于已设置并在建模中考虑热界面材料，则忽略接触热阻和接触电阻.

3）冷端温度固定为-10 ℃.

4）冷端和热端的表面对流传热系数设定为

150 W/（m²·K）， 环境温度设定为10 ℃.

5）在终端和接地中电流守恒，上、下陶瓷片都设定为电绝缘.

固定约束和机械压装是热电模块的两种基本固定方式.热端和冷端可以采用不同的固定方式，表3列出了安装边界条件的不同组合.

2.4 网格划分和模型验证

将钎料层进行自由四面体网格划分，热电臂与钎料层、钎料层与电极片间的接触表面采用自由三角形网格划分.为了使有限元仿真结果更为精确，保证计算网格有更好的质量，在交界处进行网格加密.以热电模块的剖面温度分布作为判断依据，对网格进行无关性验证，图3给出了当网格数分别为13.5万、 26.1万、51.3万和89.2万时，计算得到的模块剖面温度对比，从剖面温度曲线的比较中发现这4种网格间的最大偏差小于1%，这表明该模型满足网格独立性验证.当网格数为26.1万时，进一步增加网格数不会引起温度曲线的明显变化.因此，综合考虑计算精度与计算时间后，最终网格数确定为26.1万.

表4为C-31106型热电模块参数.为了验证有限元模型的准确性，将计算结果与表4中某公司官方实验结果进行比较，其中热端温度分别固定为25 ℃和50 ℃.

在相同条件下，改变工作电流I的大小，获得热电模块的不同制冷量Q_c，如图4所示.由图4可知，当该计算模型的热端温度T_h为25 ℃，最大电流I_max约为2.8 A时，对应的模块制冷量Q_c为5.802 W；当该模型T_h为50 ℃，I_max约为2.8 A时，对应的模块制冷量Q_c为6.308 W. 对比表4可知，当T_h分别为25 ℃和50 ℃时，其计算结果与实验结果相比，误差分别为1.66%和1.44%，说明有限元模型能较好地模拟该型热电模块实际运行时的工况.

3 计算结果分析

3.1 工作电流对模块温度的影响

为了研究热电模块在不同工作电流下的传热过程，图5给出了不同工作电流下，热电模块从热端到冷端的剖面温度（沿z轴）对比图.

由图5可见，随着工作电流的增加，模块热端温度也增大，当工作电流I从2.0 A增加到2.6 A时，热端温度T_h从295 K增加到316 K，热电臂底部（热端附近）出现温度极值点，沿z轴正方向，温度先升高后降低，可以发现当温度减小时斜率绝对值呈增大的趋势，在2.65 mm附近斜率绝对值达到最大.在z轴剖面2.50~2.65 mm范围内对比I=2.0 A和2.6 A两条曲线可以发现，电流越大温度曲线的斜率绝对值越大.说明模块产生的焦耳热在热端积累，使得温度在热电臂底部小范围升高，而在冷端带走热量，使得温度沿z轴降低.一方面是由于模块中焦耳热、傅里叶热和汤姆逊热的产生和转化；另一方面，热电材料的变物性参数和计算所用的边界条件会对模块温度分布产生一定的影响.以变物性参数和电流边界条件为例开展分析，热电材料的变物性参数可参照式（17）~式（22）， 这些参数随温度变化而变化，当它们发生变化时，反过来会影响温度的变化；工作电流的增加使得剖面温度曲线的变化较为明显.因此，工作电流越大，剖面温度曲线变化越明显，越靠近模块冷端相应位置的温度梯度也越大.

3.2 热电臂高度对模块的影响

3.2.1 热电臂高度对模块温度的影响

本节研究热电模块在不同热电臂高度下的传热过程. 在相同电流条件下，改变热电臂高度，获得热电模块从下到上各结构的剖面温度（沿z轴）对比图，如图6所示. 为了便于讨论，用模块高度H来表示热电臂高度h的变化情况.

由图6可知，随着模块高度H的增加，模块热端温度明显增大，模块冷端温度略微增大，当H从3.15 mm增加到3.95 mm时，其热端温度T_h从291 K增加到321 K，热电臂底部（热端附近）出现温度极值点，沿z轴正方向，温度先升高后降低，可以发现当温度减小时斜率绝对值逐渐增大. 对比H=3.15 mm和H=3.95 mm两条曲线可以发现，随着热电臂高度的增加，整个模块剖面温度明显升高. 这主要由热电材料的电阻值和热导值的改变所导致. 一方面，由电阻值的定义式R=ρh/A可知，假定热电臂截面面积A和材料电阻率ρ保持不变，增大热电臂高度h显然会导致其电阻增加，电阻增加又会导致焦耳热增多，恶化传热过程，使得模块整体温度升高. 另一方面，由热导值的定义式K=λA/h可知，假定A和材料热导率λ保持不变，增大h会直接导致其热导减小，从而削弱热电模块的制冷效果. 因此，热电臂高度越小，对应位置的剖面温度分布相对越均匀，温度梯度越小.

3.2.2 热电臂高度对模块应力的影响

针对不同热电臂高度对模块整体热应力场的影响进行仿真计算. 由于热电模块中的每一个热电偶都具有相同热特性^［37］，并且热电模块失效多是由于热电偶上的热应力集中及其应力最大值造成的裂纹扩展、断裂、解理性破坏^［38］. 因此，本文以模块中某一热电偶内热电臂中最大von Mises应力表征其模块内部整体热应力水平.

图7为模块高度H对冷端最大von Mises应力σ_c，max的影响.由图7可见，在4组工作电流下，当H从2.95 mm增加到3.95 mm时，热电臂冷端最大von Mises应力σ_c，max都呈现出先减小后增大的趋势，且都存在一个应力极小值. 当I为2.0~2.6 A时，其热应力极小值分别为8.98 MPa、8.68 MPa、9.15 MPa、9.77 MPa.

图8为模块高度H对热端最大von Mises应力σ_h，max的影响. 由图8可见，在不同工作电流下，热电臂热端最大von Mises应力σ_h，max随H的增大呈现先减小后迅速增大的趋势，且都存在一个应力极小值. 当I为2.0~2.6 A时，其热应力极小值分别为9.68 MPa、3.15 MPa、3.78 MPa、3.03 MPa. 对比发现，在不同工作电流下，冷、热两端都存在热应力最小值，该最小值对应一个最优热电臂高度. 因此，选择合适的热电臂高度有利于提高其热电模块的力学可靠性.

3.3 热电臂截面边长对模块的影响

3.3.1 热电臂截面边长对模块温度的影响

本节研究热电模块在不同热电臂截面边长时的传热过程. 在相同电流条件下，改变热电臂截面边长L，获得热电模块从下到上各结构位置的剖面温度（沿z轴）对比图，如图9所示. 由图9可知，随着热电臂截面边长L的增加，模块热端温度明显减小；模块冷端温度变化不大. 当L从0.89 mm增加到1.05 mm时，热端温度T_h从331 K减小到303 K，热电臂底部（热端附近）出现温度极值点，沿z轴正方向，温度先升高后降低，可以发现，当温度减小时斜率绝对值逐渐增大. 此外，对比L=0.89 mm和1.05 mm两条曲线可以发现，随着热电臂截面边长的增加，整个模块剖面温度明显降低. 由热电材料的电阻公式R=ρh/A和热导公式K=λA/h可知，假定热电臂高度h和材料电阻率ρ保持不变，增大热电臂截面边长L会直接使截面面积A增大，进而使得电阻值减小，电阻值减小又会导致焦耳热减小，有利于改善传热过程，使得模块整体温度降低. 同理，假定h和材料热导率λ保持不变，增大L会直接使截面面积A增大，进而导致热导值增加. 另外，热电臂截面边长的改变也会影响傅里叶热的大小，从而影响热电模块的温度分布，即随着L增加，在降低热阻的同时会带来傅里叶热的相对增加，从而提高热电模块的制冷性能. 因此，热电臂截面边长越大，对应位置的剖面温度分布相对越均匀，温度梯度越小.

3.3.2 热电臂截面边长对模块应力的影响

针对不同热电臂截面边长对模块整体热应力场的影响进行仿真计算. 图10为热电臂截面边长L对冷端最大von Mises应力σ_c，max的影响.由图10可见，在4组工作电流下，当L=0.89~0.99 mm时，热电臂冷端最大von Mises应力σ_c，max随L的增大而出现波动，其出现的低谷和波峰与最大von Mises应力自身的极值特性有关^［36］；当I=2.0 ~2.6 A时，在L=0.97 mm处出现了σ_c，max的极小值，分别为43.41 MPa、38.59 MPa、33.20 MPa和27.22 MPa；当L从0.99 mm增加到 1.09 mm时，σ_c，max随L的增大而明显下降，在此区间上4种电流所对应的最大von Mises 应力σ_c，max分别减小了70.13 %、76.45 %、75.77 %和74.67 %.

图11给出了热电臂截面边长L对热端最大von Mises应力σ_h，max的影响.由图11可见，在4组工作电流下，当L=0.89~1.02 mm时，热电臂热端最大von Mises应力σ_h，max随L的增大而出现先增大后上下波动；当L从1.02 mm增加到1.09 mm时，σ_h，max随L的增大呈现出逐渐下降的趋势，在此区间，4种电流所对应的最大von Mises 应力σ_h，max分别减小了27.60 %、26.30 %、24.67 %和22.60 %. 对比发现，在这一热电臂截面边长范围内（1.02~1.09 mm），热电臂截面边长越大，越有利于提高热电模块的力学可靠性.

3.4 安装方式对模块的影响

3.4.1 不同模块高度下安装方式对模块应力的影响

针对不同模块高度和安装方式对模块整体热应力场的影响进行仿真计算. 图12给出了在不同模块高度下模块冷端的最大von Mises应力σ_c，max随安装方式、机械压力P_n变化曲线. 其中，N1为上陶瓷片机械压装、下陶瓷片固定约束，N2为上陶瓷片固定约束，下陶瓷片机械压装. 由图12可见，在所有情况下，模块冷端的最大von Mises应力σ_c，max随着机械压力P_n的增大而增加. 当安装方式为N2时，σ_c，max随模块高度H的增加而增大. 以H=3.1 mm为例分析安装方式的影响，可以看到σ_c，max在采用N2的安装方式时为15.71~19.13 MPa，与采用N1时的68.09~73.46 MPa相比降低了73.95% ~76.92%.

图13给出了在不同模块高度下模块热端的最大von Mises应力σ_h，max随安装方式、机械压力P_n变化曲线. 由图13可见，当采用N1的安装方式时，σ_h，max随机械压力P_n的增大而减小.当采用N2的安装方式，且H=3.1 mm和3.3 mm时，σ_h，max随P_n的增大而增大；当采用N2的安装方式，且H=3.5 mm和3.7 mm时，σ_h，max随P_n的增大而减小. 以H=3.1 mm为例分析安装方式的影响，可以看到σ_h，max在采用N2的安装方式时为11.92~17.80 MPa，与采用N1时的49.51 ~53.65 MPa相比降低了64.05%~77.78%. 对比发现，随着机械压力的增加，采用合适的模块高度和N2的安装方式可有效提高热电模块的力学可靠性.

3.4.2 不同截面边长下安装方式对模块应力的影响

针对不同热电臂截面边长和安装方式对模块整体热应力场的影响进行仿真计算. 图14给出了在不同热电臂截面边长的条件下模块冷端的最大von Mises应力σ_c，max随安装方式、机械压力P_n变化的曲线. 由图14可见，在所有情况下模块冷端的最大von Mises应力σ_c，max随着机械压力P_n的增大而增加. 当安装方式为N1时，σ_c，max随热电臂截面边长的增加而先减小后增大. 以L=0.93 mm为例，可以看到σ_c，max在采用N2的安装方式时为25.38~32.65 MPa，与采用N1时的48.78 ~63.72 MPa相比降低了47.97%~48.76%.

图15给出了在不同热电臂截面边长的条件下模块热端的最大von Mises应力σ_h，max随安装方式、机械压力P_n变化的曲线.由图15可见，当采用N1的安装方式时，σ_h，max随机械压力P_n的增大而减小；当采用N2的安装方式时，σ_h，max随P_n的增大而增大；当采用N2的安装方式，且L=0.93 mm和0.97 mm时，σ_h，max随P_n的增大而出现更快增长.以L=0.93 mm为例，可以看到σ_h，max在采用N2的安装方式时为7.67~22.55 MPa，与采用N1时为29.53 ~22.56 MPa相比降低了74.02%~0.04%. 当P_n超过1.3 MPa后，采用N2安装方式的应力水平也开始逐渐接近并部分超过采用N1安装方式的应力水平，不过热端的热应力值仍处于较低的范围内. 对比发现，随着机械压力的增加，采用较大的热电臂截面边长和N2的安装方式可有效提高热电模块的力学可靠性.

4 结 论

建立了由31对热电偶组成的热电制冷模块多物理场模型，分析了热电臂高度、截面边长、安装方式对模块力学性能的影响，同时分析了模块内von Mises 应力的变化规律. 主要结论如下：

1）热电臂高度越小，热电臂截面边长越大，对应位置的剖面温度分布相对越均匀，温度梯度越小.

2）在对应工况下，选择合适的热电臂高度和截面边长有利于降低模块整体热应力水平，进而有效提高热电模块的力学可靠性.

3）对比发现，随着机械压力的增加，对模块热端采用机械压装的安装方式，冷端采用固定约束的方式可显著降低模块整体热应力水平，进而有效提高热电模块的力学可靠性.

参考文献

[1] 王俊强， 高利军， 李运泽. 基于TEC的空间站末端回路温控系统建模及其热力学性能分析［J］. 航空动力学报， 2019， 34（7）：1483-1492．

[2] WANG J Q， GAO L J， LI Y Z. Modeling and thermodynamic performance analysis of thermal control system based on terminal circuit with TEC in space station［J］. Journal of Aerospace Power，2019，34（7）：1483-1492．（in Chinese）

[3] 江帆，孟凡凯，陈林根，等 ．变温热源小型热电冷水机结构设计与性能分析［J］．工程热物理学报，2020，41（7）：1573-1578．

[4] JIANG F， MENG F K， CHEN L G，et al ．Structural design and performance analysis of a small thermoelectric chiller with variable temperature heat reservoirs［J］．Journal of Engineering Thermophysics，2020，41（7）：1573-1578．（in Chinese）

[5] 陈林根， 孟凡凯， 戈延林， 等 ．半导体热电装置的热力学研究进展［J］. 机械工程学报，2013，49（24）： 144-154．

[6] CHEN L G， MENG F K， GE Y L，et al. Progress in thermodynamic studies for semiconductor thermoelectric devices［J］．Journal of Mechanical Engineering，2013，49（24）：144-154．（in Chinese）

[7] GONG T R， GAO L， WU Y J， et al. Numerical simulation on a compact thermoelectric cooler for the optimized design［J］．Applied Thermal Engineering， 2019， 146： 815-825．

[8] WANG C Y， YANG X， SHEN Y N， et al ．Numerical investigation on cooling performance of multilayer pyramid thermoelectric module［J］. International Communications in Heat and Mass Transfer，2023，143：106738．

[9] LI S， LIU J L， DING L，et al ．Active thermal management of high-power LED through chip on thermoelectric cooler［J］．IEEE Transactions on Electron Devices，2021，68（4）： 1753-1756．

[10] CHU W S， LÜ J H， HONG N ．Finite element analysis of Peltier effect based thermoelectric energy conservation cooling system for microprocessors［J］. IOP Conference Series：Earth and Environmental Science，2021，804（3）：032056．

[11] 王镇锐，张兴斌，温世喆，等 ．结合TEC的泵驱两相温控系统的空间应用［J］．宇航学报，2018，39（10）：1176-1184．

[12] WANG Z R， ZHANG X B， WEN S Z，et al ．Space applications of pumped two-phase temperature control system combined with TEC［J］．Journal of Astronautics，2018，39（10）：1176-1184．（in Chinese）

[13] 宋雪冬，周琦，李猛，等 ．基于FPGA的星敏感器探测器制冷控制系统［J］．飞控与探测，2022，5（5）： 47-51．

[14] SONG X D， ZHOU Q， LI M，et al ．Refrigeration control system of star sensor detector based on FPGA［J］. Flight Control & Detection，2022，5（5）：47-51．（in Chinese）

[15] MENG F K， CHEN L G， FENG Y L，et al ．Thermoelectric generator for industrial gas phase waste heat recovery［J］．Energy，2017，135：83-90．

[16] ZHAO L J， YANG J， ZOU Y H，et al ．Tuning Ag content to achieve high thermoelectric properties of Bi-doped p-type Cu₃SbSe₄-based materials［J］．Journal of Alloys and Compounds，2021，872：159659．

[17] 徐庆， 赵琨鹏， 魏天然， 等. 热电材料的研究现状与未来展 望［J］. 硅酸盐学报， 2021， 49（7）： 1296-1305．

[18] XU Q， ZHAO K P， WEI T R， et al ．Development and prospects of thermoelectric materials［J］．Journal of the Chinese Ceramic Society，2021，49（7）：1296-1305．（in Chinese）

[19] CHEN L G， MENG F K， XIE Z H，et al ．Thermodynamic modeling and analysis of an air-cooled small space thermoelectric cooler［J］. The European Physical Journal Plus，2020，135（1）：80．

[20] BALDRY M， TIMCHENKO V， MENICTAS C ．Optimal design of a natural convection heat sink for small thermoelectric cooling modules［J］．Applied Thermal Engineering，2019，160：114062．

[21] ZHU L， YU J L ．Optimization of heat sink of thermoelectric cooler using entropy generation analysis［J］．International Journal of Thermal Sciences，2017，118：168-175．

[22] MENG F， CHEN L， SUN F ．Performance prediction and irreversibility analysis of a thermoelectric refrigerator with finned heat exchanger［J］．Acta Physica Polonica A，2011，120（3）：397-406．

[23] 王子成，蔡兰兰，高鹏，等 ．热电制冷强化风冷散热模块的工作特性分析［J］．制冷学报，2020，41（2）：48-55．

[24] WANG Z C， CAI L L， GAO P，et al ．Operating characteristics analysis of thermoelectric cooler enhanced air cooling module［J］．Journal of Refrigeration，2020，41（2）：48-55．（in Chinese）

[25] 孟凡凯，徐辰欣，孙悦桐 ．密闭空间串联式两级热电制冷器瞬态特性分析［J］．湖南大学学报（自然科学版），2023，50（6）：99-107．

[26] MENG F K， XU C X， SUN Y T ．Analysis on transient characteristics of series two-stage thermoelectric cooler in confined space［J］．Journal of Hunan University （Natural Sciences），2023，50（6）：99-107．（in Chinese）

[27] ANG E Y M， NG P S， SOH C B，et al ．Multi-stage thermoelectric coolers for cooling wearables［J］．Thermal Science and Engineering Progress，2022，36：101511．

[28] 孙悦桐，孟凡凯，周林，等 ．空冷式多级热电制冷器性能综合分析［J］．上海交通大学学报，2024，58（3）：371-381．

[29] SUN Y T， MENG F K， ZHOU L，et al ．Comprehensive analysis of performance of air cooled multistage thermoelectric cooler［J］．Journal of Shanghai Jiao Tong University，2024，58（3）：371-381．（in Chinese）

[30] 毛东雪，汪文帅， 谢军 ．远场载荷下含椭圆孔热电薄板的热电力耦合解析解［J］．工程力学，2023，40（4）：21-34．

[31] MAO D X， WANG W S， XIE J ．Thermoelectric coupled analytical solution of thermoelectric plate with elliptical hole subjected to far field load［J］．Engineering Mechanics， 2023， 40（4）： 21-34．（in Chinese）

[32] 周岷峰 ．热电有限元模型应力分析研究现状及方法介绍［J］．中国科学（技术科学），2013，43（4）：397-406．

[33] ZHOU M F ．Present situation and method introduction of stress analysis of thermoelectric finite element model［J］. Scientia Sinica （Technologica），2013，43（4）：397-406．（in Chinese）

[34] 刘冠宇，王红梅，黄荣森，等 ．热电制冷器的结构优化研究［J］．制冷技术，2021，41（1）：55-61．

[35] LIU G Y， WANG H M， HUANG R S，et al ．Research on structural optimization of thermoelectric cooler［J］. Chinese Journal of Refrigeration Technology，2021，41（1）：55-61．（in Chinese）

[36] BARAKO M T， PARK W， MARCONNET A M，et al ．Thermal cycling，mechanical degradation，and the effective figure of merit of a thermoelectric module［J］．Journal of Electronic Materials，2013，42（3）：372-381．

[37] 周林， 孟凡凯， 陈华伟， 等. 钎料层对热电模块力学性能的影响分析［J/OL］. 工程力学. （2023-02-11）［2023-08-27］.

[38] ZHOU L， MENG F K， CHEN H W， et al. Analysis on influence of solder layers on mechanical properties of a thermoelectric module［J/OL］. Engineering Mechanics. （2023-02-11）［2023-08-27］. in Chinese）

[39] ARANGUREN P， ASTRAIN D， RODRÍGUEZ A，et al ．Net thermoelectric power generation improvement through heat transfer optimization［J］．Applied Thermal Engineering，2017，120：496-505．

[40] TAN H B， FU H， YU J L ．Evaluating optimal cooling temperature of a single-stage thermoelectric cooler using thermodynamic second law［J］．Applied Thermal Engineering，2017，123：845-851．

[41] 李俊伟 ．方钴矿基热电器件的应力仿真研究［D］．上海：中国科学院大学，2021．

[42] LI J W. Stress simulation of skutterudite-based thermoelectric devices［D］. Shanghai：University of Chinese Academy of Sciences，2021．（in Chinese）

[43] ZHOU L， MENG F K， SUN Y T ．Numerical study on infrared detectors cooling by multi-stage thermoelectric cooler combined with microchannel heat sink［J］．Applied Thermal Engineering，2024，236：121788．

[44] CONFALONIERI C， PERRIN M， GARIBOLDI E ．Combined powder metallurgy routes to improve thermal and mechanical response of Al-Sn composite phase change materials［J］．Transactions of Nonferrous Metals Society of China， 2020， 30（12）： 3226-3239．

[45] OOTAO Y， TANIGAWA Y ．Three-dimensional solution for transient thermal stresses of functionally graded rectangular plate due to nonuniform heat supply［J］．International Journal of Mechanical Sciences，2005，47（11）：1769-1788．

[46] 张涛然， 晁晓洁， 郭丽红. 材料力学［M］. 重庆： 重庆大学出版社， 2018： 148-171.

[47] ZHANG T R， CHAO X J， GUO L H. Mechanics of materials［M］. Chongqing： Chongqing University Press， 2018： 148-171. （in Chinese）

[48] 张子豪 ．纳米结构Bi_0.5Sb_1.5Te₃热电材料的连接工艺及机理研究［D］．哈尔滨：哈尔滨工业大学，2017．

[49] ZHANG Z H ．Processing and mechanism of joining Bi_0.5Sb_1.5Te₃ thermoelectric material with nanostructure［D］．Harbin：Harbin Institute of Technology，2017．（in Chinese）

[50] BAI N， CHEN X， GAO H ．Simulation of uniaxial tensile properties for lead-free solders with modified Anand model［J］．Materials & Design，2009，30（1）：122-128．

[51] SON J S， CHOI M K， HAN M K，et al ．N-Type nanostructured thermoelectric materials prepared from chemically synthesized ultrathin Bi₂Te₃ nanoplates［J］．Nano Letters，2012，12（2）：640-647．

[52] JIANG F， MENG F K， CHEN L G，et al ．Thermodynamic analysis and experimental research of water-cooled small space thermoelectric air-conditioner［J］．Journal of Thermal Science，2022，31（2）：390-406．

[53] 范世发 ．几类热电器件的热-电-力学性能分析及结构优化设计［D］．兰州： 兰州大学，2021．

[54] FAN S F. Thermal-electro-mechanical performances analysis and structural optimization design of several types of thermoelectric devices［D］. Lanzhou：Lanzhou University，2021．（in Chinese）

[55] GAO J L， DU Q G， ZHANG X D，et al ．Thermal stress analysis and structure parameter selection for a Bi₂Te₃-based thermoelectric module［J］．Journal of Electronic Materials，2011，40（5）：884-888．

基金资助

国家自然科学基金资助项目（11974429）; 海军工程大学自主研发计划项目（2022502210）